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ABSTRACT 

We present a new approach to investigate the content and spatial distribution of dust in struc- 
turally unresolved star-forming galaxies from the observed dependence of integrated spectral 
properties on galaxy inclination. Motivated by the observation that different stellar populations 
reside in different spatial components of nearby star-forming galaxies, we develop an inno- 
vative combination of generic models of radiative transfer in dusty media with a prescription 
for the spectral evolution of galaxies, via the association of different geometric components of 
galaxies with stars in different age ranges. We start by showing that a wide range of radiative 
transfer models all predict a quasi-universal relation between slope of the attenuation curve 
at any wavelength, from the ultraviolet to the near-infrared, and V-band attenuation optical 
depth in the diffuse interstellar medium (ISM), at all galaxy inclinations. This relation pre- 
dicts steeper (shallower) dust attenuation curves than both the Calzetti and Milky-Way curves 
at small (large) attenuation optical depths, which implies that geometry and orientation effects 
have a stronger influence on the shape of the attenuation curve than changes in the optical 
properties of dust grains. We use our new, combined radiative transfer and spectral evolution 
model to interpret the observed dependence of the Ha/H/3 ratio and ugrizYJH attenuation 
curve on inclination in a sample of about 23 000 nearby star-forming galaxies, which we 
correct for systematic biases by developing a general method based on importance sampling. 
From the exploration of the model parameter space by means of a Bayesian Markov Chain 
Monte Carlo technique, we measure the central face-on B-band optical depth of this sample 
to be Tg± ~ 1.8 + 0.2 (corresponding to an angle-average (f^ M )e » 0.3). We also quantify 
the enhanced optical depth towards newly formed stars in their birth clouds, finding this to 
be significantly larger in galaxies with bulges than in disc-dominated galaxies, while tb± is 
roughly similar in both cases. This can arise if, for example, galaxies with significant bulges 
have higher central star formation efficiencies than their disc-dominated counterparts at fixed 
specific star formation rate, and dustier stellar birth clouds because of the higher metallicity. 
We find that over 80 percent of the attenuation in galaxies in our sample is characteristic of 
that affecting thin-disc stars in radiative transfer models. The median unattenuated V-band 
luminosity ratio of thick-disc to thin-disc stars is 0.1-0.2, in good agreement with the results 
from spatially resolved studies of nearby edge-on disc galaxies. Finally, we show that neglect- 
ing the effect of geometry and orientation on attenuation can severely bias the interpretation 
of galaxy spectral energy distributions, as the impact on broadband colours can reach up to 
0.3-0.4mag at optical wavelengths and 0.1 mag at near-infrared ones. This paper also con- 
tains an original application of Gaussian Random Processes to extend the wavelength range 
of dust attenuation curves. 

Key words: galaxies: general, galaxies: ISM, dust, extinction, radiative transfer 



1 INTRODUCTION 

A small but significant mass fraction of the interstellar medium 
(ISM) in galaxies is in the form of dust grains, which abso rb and 
scatte r the light emitted by stars at all wavelengths (e.g. Spitzer 
* e-mail:chevalla@iap.fr 1 1978b . Dust grains are produced during late phases of stellar evolu- 
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tion - mainly supernova ex plosions and winds of asymptotic-giant- 
branc h (AGB) stars (e.g. |Pweklll998l ; iHofnerl 120091 ; Icherchnefj 
2010) - and are destroyed by energetic photons and shocks (e.g. 
Jonesl 120041 : iTielenj I200II ; I Jones & NutH 1201 ll) . This makes the 
amount, composition and spatial distribution of dust in a galaxy 
depend in a complex way on the star formation history and geom- 
etry, and hence, the history of hierarchical merging, gas infall and 
outflow. All these factors influence the wavelength dependence of 
the attenuation of starlight by dust, which we must understand to 
retrieve star formation and chemical enrichment histories from ob- 
served spectral energy distributions of galaxies. 

Much information has been gathered on the shape of the dust 
attenuation curve in galaxies. In the Milky Way, the Large and 
the Small Magellanic Clouds (LMC and SMC), comparisons of 
observed spectral energy distributions of stars of fixed spectral 
type and luminosity class along different lines of sight have pro- 
vided esti mates of the extinction*! curve along individual lines of 
sight (e.g.lSavage & Mathisll979flPrevot et aljl9 84; Bouch et et al] 
ll985tlClavton & Martinll 19851; ICardelli et al J 1 1 98 8L U989I) . By ex- 
tending this approach to the analysis of the integrated light in a 
sample of nearby starburst galaxies, for which the unattenuated 
spectral ener gy distribution is k nown to be dominated by blue 
massive stars, ICalzetti et all i ll 994 see also lCarzettfeOOlh derived 
constraints on the mean attenuation curve in these galaxies. This 
is significantly greyer than the Milky Way, LMC and SMC ex- 
tinction cur ves at ultraviol et wavelengths. Theoretically, the grey- 
ness of the lCalzettil d200lh attenuation curve can be accounted for 
by a simple angle- averaged, two-component dust model, in which 
young stars in giant molecular clouds suffer more attenuation than 
older stars in the diffuse ISM, for standard (i.e. LM C-like) op- 
tical properties of dust grains dCharlot&Falll boOO). Models of 
this type are useful because they can be implement ed easily in 
any spectral synthesis analy s is of galaxies (e.g. 



Brinchmann et alj 



|2004 Ida Cunha et al. 2008; |PacificietalJ [2012). Recently, how- 
ever, the identification of systematic changes of the dust attenua- 
tion curve as a function of both inclination and star formation his- 
tory in nearb y galaxies has show n the limitation of angle-averaged 
dust models dWild et alj|201 lbl) . Investigating the physical origin 
of these changes requires more sophisticated models to compute 
the transfer o f photons through the dusty ISM of spatially resolved 
galaxies (e .g. Silva et alJl998l ; lTuffs et alj2004|Pierini et aljiooj ; 



Ijonsson et al. l2010h . In practice, such models are computationally 
expensive and cannot be easily applied to the detailed analysis of 
large samples of galaxies. 

In this paper, we present a new approach to account in a simple 
yet physically consistent way for the dependence of dust attenua- 
tion on orientation in spectral analyses of structurally unresolved 
galaxies. This approach is based on an original exploitation of ex- 
isting sophisticated models of radiative transfer in dust y media. We 
consi der four types of models rel ying on Monte-Ca rlo dPierini et alj 
2004 hereafter P04). analytic dTuffs et al.1 12004 hereafter T04; 
Silva et al] 1 19981, hereafter S9 8), smoothed particle hydrodynam- 
ics (SPH. Ijonsson et al] I2010L hereafter J 10) calculations of the 
scattering and absorption of photons through the diffuse ISM of 
galaxies with multiple-disc and bulge stellar components. We show 



1 We follow the standard nomenclature that 'extinction' refers to the ab- 
sorption of photons along, and their scattering out of, a line of sight, while 
'attenuation' (or 'effective absorption') refers to the combined effects of ab- 
sorption and scattering in and ou t of the line of sight caused by both local 
and global geometric effects (e.g. Calzetti 2001). 



that these different types of models all predict a quasi-universal re- 
lation between slope of the attenuation curve at any wavelength 
from the ultraviolet to the near infrared and V-band attenuation 
optical depth of the dust in the diffuse ISM, at all galaxy incli- 
nations. We further show that appropriate choices of the dust op- 
tical depths and relative luminosities of the various disc and bulge 
components in the flexible T04 model can reproduce the depen- 
dence of attenuation on galaxy inclination predicted by the P04 and 
J10 models. On these stable grounds, we combine the T04 model 
with a semi-an alytic post-treatment of the Millennium cosmo logi- 
cal simulation (Sprin gel et al.l2005l;lDe Lucia & Blaizotl2007h . the 
latest version of the iBruzual & Charlotl <l2003r) stellar population 
synthes is code and a prescrip tion for dust attenuation in stellar birth 
clouds dCharlot & FaTfo oO). to investigate the physic al origin of 
the obs erved changes of the dust attenuation curve in the lWild et alj 
d201 lbl) sample of nearby galaxies. We achieve this by relating the 
geometric components of galaxies to stellar age ranges in spatially 
unresolved models of spectral evolution. In practice, we build a 
comprehensive library of model spectral energy distributions en- 
compassing wide ranges of stellar and dust parameters for galaxies 
with multiple-disc and bulge components. Then, we use this library 
to constrain the dust content and stellar age ran ges in the differen t 
geometric components of the galaxies in the IWild et aT] d201 lbh 
sample by performing a Bayesian analysis of the observed Ha/HfS 
ratio and broadband (optical and near-infrared) spectral energy dis- 
tribution. The results of this analysis provide a calibration of the 
dependence of attenuation on both orientation and star formation 
history in structurally unresolved galaxies. We show that account- 
ing for this effect has important implications for the interpretation 
of galaxy spectral energy distributions. 

In Section 2 below, we briefly recall the formalism to describe 
the attenuation of starlight by dust in galaxie s and compare the ra - 
diative transfer models of iTuffs et alj d2004) . IPierini et al] d2004) . 

t et al.1 dll 



Ijonsson et aT] d2010h and lSilva et al.1 dl998l) . We show that the pre 



dictions of these different models are consistent with one another 
and yield a quasi-universal relation between slope of the attenua- 
tion curve at any wavelength and attenuation optical depth of the 
dust. In Section 3, we combine the T04 model with models of the 
spectral evolution of galaxies in a cosmological context to build a 
comprehensive library of synthetic spectral energy distributions of 
dusty galaxies. We use this l ibrary to int e rpret th e observed line and 
continuum emission of the IWild et aT] d201 lbl) sample of nearby 
star-forming galaxies in terms of constraints on the content and dis- 
tribution of the dust between the different stellar components. We 
discuss the implications of these results for the analysis of struc- 
turally unresolved galaxies in Section 4, where we also summarise 
our conclusions. 



2 GENERAL PROPERTIES OF DUST ATTENUATION 
MODELS 

The intensity I A {S) emerging at wavelength A in a direction 9 from 
the normal to the equatorial plane of a galaxy can always be re- 
lated to the intensity if produced by stars, assumed isotropic, by an 
expression of the form 



I A (8) = I°expl-U0)] , 



(1) 



where f A (d) is the 'attenuation' (or 'effective absorption') optical 
depth of the dust affecting those photons emitted in all directions by 
all stars in the galaxy that emerge in the direction 6. Expression [TJ 
implicitly assumes azimuthal symmetry. In practice, the quantity 
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f A {8) depends on the geometry of the galaxy and the optical prop- 
erties and spatial distribution of dust relative to the stars. This can 
make the intensity 1 A {8) either smaller or larger than 7°, depending 
on the influence of scattering on the paths of photons before these 
escape from the galaxy. 

In practice, the attenuation optical depth f A (8) incorporates 
contributions from the different geometric components of the 
galaxy (i.e. bulge, thin and thick discs), each influencing in a dif- 
ferent way the dependence of f A on 8. We write 



f,(e) = -ln^^ 



exp[-f ,,,«?)]. 



(2) 



where f Aj {8) is the contribution by the i th component to the absorp- 
tion optical depth of the dust affecting photons emerging at wave- 
length A in the direction 8, and 



(3) 



is the contribution by this component to the intensity produced by 



stars at wavelength A (with £,■/, 
length A can be expressed as 



/ ). The intensity at wave- 



Il(t) 



f 

Jo 



dt'tfr(t-t')S A [t',Z(t-t')], 



(4) 



where i[/(t-f) is the star formation rate at time t — f and S A [t',Z(t - 
f')] is the intensity produced per unit wavelength and per unit mass 
by a stellar population of age f and metallicity Z(t - t'). 

For some purposes, it is useful to compute the mean of f A (9) 
over all solid angles. In the case of azimuthal symmetry (dQ = 
2n rfcos 8), and if we assume that the radiation is symmetric with 
respect to the equatorial plane of the galaxy, this is given by 



(We : 



f 

Jo 



dcos 6 f A (8) . 



(5) 



2.1 Radiative Transfer models 

The attenuation optical depth f A (8) can be computed analytically 
for simple geometries, such as a uniform screen of dus t in front 
of the stars and a uniform mixture of stars and dust (e.g. lYip et al] 
|2010|) . The study of more complex geometries requires a numeri- 
cal approach to follow the interactions between photons and dust 
grains in the interstellar medium. We focus here on f our recent, 
popula r models of radiativ e transfer in dusty m edia: 
d2004). IPierini et al.1 | |2004|) . Ijonsson et all feOld) and 



Tuffs et al 



Silva et al 



(1998). Fig. [T] shows that the dependence on wavelength of the 
dust-absorption coefficient, single-scattering albedo and asymme- 
try parameter of the scattering phase function adopted in these dif- 
ferent models for Milky-Way type dust are roughly consistent with 
one another, except for the fact that S98 assume isotropic scatter- 
ing at all wavelengths (i.e. g A = 0). The small apparent discrepan- 
cies between the dust optical properties of the P04, T04 and J 10 
models have a negligible impact on our analysis. In fact, we show 
in Sections |2.2| and [231 below that the attenuation curves obtained 
from these radiative transfer models have similar global properties 
at fixed geometry and inclination, while the adoption of isotropic 
scattering in the S98 model introduces significant differences in the 
shape of the attenuation curve at low attenuation optical depths. In 
the next paragraphs, we briefly review the T04, P04, J10 and S98 
models and establish their common behaviours and mutual coher- 
ence. We focus here on the predictions of these models for atten- 
uation by dust in the ambient (i.e. diffuse) ISM, which we refer to 



3 0.5 




Figure 1. Dust-absorption coefficient in units of the V-band value, k a /kv, 
asymmetry parameter of the scattering phase function, g A , and single- 
scattering albedo, co A , plotted against wavelength A, as adopted for Milky - 
Way type dust in th e P04 (short-dashed line; from | Witt & Gordonll2000h . 
T04 (solid line; fromlLaor & Draine|[ l993; Math is et alj 19771) . J10 (dotted 
line; from lWeingartner & DraindbOOlt famine & Lill2007l) and S98 (long- 
dashed l ine; fromlDraine & Led 1984 lLaor & Drainelll993l with modifica- 
tions bv lSilva et al.ll998l) models. 



Table 1. Summary of the parameters defining the geometric configuration 
of the different radiative transfer models described in Section 12.11 



Model 


Type 


scale length (kpc ) 


scale height (kpc ) 






stars (ft,) 


dust (hj) 


stars (z„) 


dust (id) 


P04 


Bulge 
Disc 


1.0 

3.0 


3.0 


0.06-0.37 


0.11 


T04 a 


Bulge 
Thin Disc 
Thick Disc 


0.69 
3.0 
4.2 


3.0 
4.2 


0.048 
0.22 


0.048 
0.14 


J10 b 


G-series 
Sbc-series 


1.1-2.8 
4.0-7.0 


3 h, 
3 ft, 


0.125 -0.2 ft, 
0.125-0.2 ft, 




S98 


S-series 
N-series 


5 

5-8 


h, 
h, 


0.5 
0.4-0.1 


Z. 

z. 



1 T04 express all parameters in units of the thin-disc scale length, while P04 
and J10 adopt physical units. For the sake of comparison, we re-normalise 
here the T04 parameters to the scale length of the P04 disc. 

' J10 do not provide any measure of the dust disc scale height, which is very 
sensitive to the resolution of the adopted SPH simulation for each galaxy. 



as t\ (8). The characteristics of these models are summarised in 
Table □ 



2.1.1 The P04 model \Pierini et alXlOOA) 

The P04 model rel ies on the DIRTY radiative transfer code of 
iGordon et alj d200ll) . This follows the multiple scattering and ab- 
sorption of photons through an analytic distribution of dust and 
stars by means of a Monte Carlo approach. P04 use DIRTY to com- 
pute the attenuation optical depth f' SM (#) for galaxies consisting 
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Stellar disc of wavelength- 
dependent scale height Bulge 




Dust disc 



Figure 2. Geometric components of the radiative transfer model of P04: 
a stellar disc with a double exponentially declining light profile and 
wavelength-dependent scale height (blue-to-red rectangle); a stellar bulge 
with an exponentially declining radial light profile (red circle); both pene- 
trated by a dust disc with a double exponentially declining density profile 
(white diagonal hatching). See Section l2. 1.1 l and Table[7Jfor detail. 



Thick stellar disc 



Bulge Thin stellar disc 
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Thick dust disc 



Thin dust disc 



Figure 3. Geometric components of the radiative transfer model of T04: a 
'thin' and a 'thick' stellar disc with a double exponentially declining light 
profile (cyan and blue rectangles, respectively); a stellar bulge with a de 
Vaucouleurs exp(-r'' 4 ) declining radial light profile (red circle); each pen- 
etrated by a 'thin' and a 'thick' dust discs with a double exponentially de- 
clining density profile (red diagonal hatching). The aspect ratios of the dif- 
ferent components reproduce those adopted in the model. See Section l2.1.2l 
and Table[T]for detail. 



of a stellar bulge and disc pervaded by a dust disc, as depicted in 

Fig.m 

The 'stellar bulge' is modelled as a sphere with an exponen- 
tially declining radial light profile of scale length 1 kpc truncated at 
4kpc. The 'stellar disc' has a double exponentially declining light 
profile of the form 

p(rJZ|)<xexpf-~-^j i (6 ) 

where r and \Z\ indicate the radial and polar coordinates and h 
and z the scale length and scale height of the disc, respectively. 
P04 adopt a fixed scale length for the stellar disc, h, = 3 kpc. To 
account for the fact that young (blue) stars tend to have a lower 
scale height than older (redder) stars in observed spir al galaxies 
(e.g. lYoachim & Dalcanton|[20oi : lYoachim et ai]l2012l) . P04 take 
the scale height of the stellar disc to increase with wavelength, from 
z, = 60 pc at A = 0.1 ixm to 370 pc at 3 ixm. The 'dust disc' also 
has a double exponentially declining density profile, with the same 
scale length as the stellar disc (hj = 3 kpc) and a fixed scale height 
Zd = HOpc, which equals that of the stellar disc at wavelength 
A = 1850 A. The stellar and dust discs both extend to a maxi- 
mum radius of 12 kpc and maximum height of 2.1 kpc. P04 com- 
pute f' SM (#) for two types of dust distribution: homogeneous and 
clumpy0 We consider only the clumpy distribution here because 
it better mimics the patchy and filamentary structure of the ISM 
(s ee, e.g. the recent observations with the Hersch el/PACS camera 
in Mookeri ea et alj20lTI ; lArzoumanian etal 1l201ll) . 

P04 conveniently parametrize their results in terms of the cen- 
tral face-on V-band optical depth of the dust from the equatorial 
plane to the surface of a galaxy, which we denote by T^f. They 
consider values of this parameter in the range 0.25 ^ T™'f ^ 8.00 
and compute f' SM (#) at angles in the range 0° ^ 9 ^ 90° and 
wavelengths in the range 0.1 A ^ 3.0 urn. We note that, even 
for T™f = 8.0, the average over angles <Ty M } e amounts to only 1.6 
because of the exponentially declining density profile of the dust 
disc. 

2.1.2 TheT04 model iTuffs et al\2004) 

The T04 model relie s on the radiative transfer code of 
iKvlafis & Bahcalj ll987T) . As in the case of the P04 models, this 

2 They describe the clumpy medium by associating in a stochastic way 
each resolution element (corresponding to a cubic cell of L = 44 pc on the 
side) to either a high or a low density state. The high-to-low density contrast 
is set to 100 and the filling factor of the high density clumps to 0.15. 



code includes multiple anisotropic scattering and can be applied 
to analytic distributions of stars and dust. Unlike the Monte Carlo 
approach of P04, the code adopted by T04 is based on an approxi- 
mate solution to the radia tive transfer equation, wh ich can be easily 
computed (for detail, see IKvlafis & Bahcail 19871) . 

T04 describe galaxies as the sum of three stellar components 
(a thin and a thick disc, and a bulge) and two dust components (as- 
sociated with the stellar discs), as illustrated by Fig. [3] The two 
stellar discs, along with the associated dust discs, are described as 
double exponentially declining light and density profiles, respec- 
tively (equation [6j. T04 normalise all the scale parameters of a 
galaxy to the fi-band scale length of the thin stellar disc, which is 
by construction /?*'" = 1.0. The 'thin stellar disc' has a scale height 
zf m = 0.016 and is pervaded by a 'thin dust disc' of same scale 
height zj m and scale length ft, . The 'thick stellar disc' has larger 
scale length and scale height than the thin disc (ft* ick = 1.406, 
z thick _ 0.074 ss 4.6 4 hm ) and it is penetrated by a dust disc with 
the same scale length and smaller scale height than the thick stellar 
disc (z?"* ~ 0.65z* lck ). The 'stellar bulge' has a de Vaucouleurs 
exp(r 1/4 ) light profile, at variance with the exponentially declining 
light profile in the P04 model. We stress that each stellar component 
in the T04 model (bulge, thick and thin disc) is attenuated by both 
the thin and the thick dust discs, but with different implications for 
f^ M (6) because of the different relative distributions of stars and 
dust. 

T04 parametrize the dust content of a galaxy by means of 
the central face-on 6-band optical depth, which includes the con- 
tributions from the two dust discs, t b± = r*} ck + t*™. They fix 
T B± k / T B'± = °- 387 and explore the range 0.1 «C r s± sC 8.0. T04 
compute tf M (9) at angles in the range 0° ^ 6 ^ 90° and wave- 
lengths in the range 0.09 ^ A ^ 2.2 pin for the thin disc, but only 
in the range 0.45 ^ A ^ 2.2 urn for the thick disc and bulge. While 
stars in the thick disc and bulge are not expected to dominate the 
radiation at wavelengths 0.09 ^ A ^ 0.45 urn, the lack of informa- 
tion on t l f M (6) in this range limits the exploration of the attenua- 
tion of potentially interesting populations of evolved stars in these 
two components. To overcome at least partially this limitation, we 
can try to guess the attenuation at A ^ 0.45 p.m based on that at 
A ^ 0.45 urn for the thick stellar disc and bulge, using the infor- 
mation available at all wavelengths for the thin stellar disc0 This 



3 Unfortunately, it was not possible for the authors of T04 to compute 
flS M (0) j n me m ; ssm g wavelength domain for our purpose. 
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is possible because, as noted above, the same dust discs determine 
the attenuation of all three stellar components. 

To best-guess the shape of the attenuation curve fJ, SM (0) at 
wavelengths A ^ 0.45 \im from that at A ^ 0.45 \im for the 
thick stellar disc and bulge, we do not choose extrapolation meth- 
ods invoking parametric functions (e.g. polynomials, power laws). 
Instead, we use a 'Gaussian random process' to learn the shape of 
the attenuation curves at A ^ 0.45 \im from that at A ^ 0.45 \im 
based on the calculations available for the thin-disc attenuation (see 
Appendix [A] for detail). This approach is valid under two main as- 
sumptions: firstly, that the attenuation curves red- ward and blue- 
ward of 0.45 |J.m be strongly enough correlated that the knowl- 
edge of one allows the prediction of the other; the finding by 
Carde lTTet al.l d 19881) that the ultraviolet extinction law correlates 
well with the optical/near-infrared one for stars in various Galac- 
tic environments supports this hypothesis. Secondly, that our learn- 
ing sample (the thin-disc attenuation curves for various t b ± and 9) 
well cover the parameter space sampled by the incomplete atten- 
uation curves (of the thick disc and bulge). In our case, the thin- 
disc calculations of T04 do encompass the full range of optical 
depths sampled by the thick-disc and bulge calculations. In Ap- 
pendix[A] we describe and thoroughly test our algorithm to extend 
the thick-disc and bulge attenuation curves computed by T04. From 
the analysis of different datasets, we show that ff M (6) at wave- 
lengths 0.25 < A ^ 0.45 Jim can be retrieved reliably from the 
attenuation curve at A ^ 0.45 (i.m, but that the recovery worsens at 
wavelengths A < 0.25 |J.m. 

We note that T04 account in a separate way for the attenua- 
tion of newly born stars in stellar birth clouds by introducing an 
analytic component of clumpy dust, which absorbs a fraction of the 
ultraviolet radiation. This component is not shown in Fig. [3] and it 
is not used as part of the T04 model in this work, since we inde- 
pendently model the attenuation of n ewly born stars following the 
prescription of ICharlot & Falll d2000l see Section |3 . 2 . 2 I belo w) , 

2. 1.3 The J 10 model Uonsson et al\20ldi) 

The J 10 model relies on the SUNRISE radiative transfer code of 
Ijonssonl fc006l) . J 10 appeal to a set of smoothed particle hydro- 
dynamics simulations to describe the input distributions of stars 
and gas in SUNRISE, which differ from the analytic profiles used 
by P04 and T04. Also, J10 adopt a Monte Carlo method to solve 
the equation of radiative transfer, as P04, but considering polychro- 
matic rays instead of monochromatic ones0 With this approach, 
single rays can be used to trace the paths of photons of differ- 
ent wavelengths, which enables high spectral resolution across the 
whole range from ultraviolet to infrared wavel engths (for detail and 
a discussion of drawbacks, see lJonsso n 200i). 

The SPH simulati ons in the J10 model were performe d with 
the GADGET code of ISpringel etail d200ll ISpringelll2005l) and 
include supernova feedback and metal enrichment dJonsson et alj 
l200d ICox et aljfeood 120081) . J 10 compute the attenuation curve 
rt SM (0) for two sets of SPH simulations: the 'G' series, with stel- 
lar and interstellar properties typical of nearby galaxies with stellar 
masses in the range 1 x 1 9 < M./M H < 6 x 10 10 in the Sloan 
Digital Sky Survey (SDSS. lYork et alfeOOOh ; and the 'Sbc' series, 
with properties typical of local late-type spiral galaxies with masses 

4 This is achieved by expressing the probability for photons of any wave- 
length to be absorbed or scattered by dust along a ray in terms of that of a 
photon of reference wavelength 0.9 p.m (see J10 for detail). 



in the range 5 x 10'° sC M,/M Q 1 x 10" (see figs 6 and 7 
of lJonsson et~aT]|201Cl) . All simulations include the self-consistent 
treatment of the evolution of dark matter, stellar and gas particles. 
Attenuation by dust is added a posteriori, assuming that a constant 
fraction of 40 percent of the mass of metals in the gas is in form 
of dust grains. In the original simulations of J 10, the attenuation of 
newly born stars in their birth cl ouds is computed using the photo- 
ionisation code MAPPINGSIII iGroves et al]|2004l , l2008h - In the 
present study, for the purpose of comparison with the P04 and T04 
models, we use a set of J 10 simulations which do not include the 
attenuation of young stars in stellar birth clouds (we are grateful 
to B. Groves and P. lonsson for kindly making these calculations 
available to us). 

A main advantage of the J 10 model over the T04 and P04 ones 
is to rely on a physically motivated spatial distribution of stars and 
dust. The drawback is that this approach is computationally expen- 
sive: the simulations of only 7 galaxies published to date do not al- 
low proper statistical analyses. Moreover, J 10 show that the results 
of their calculations depend on the resolution of the SPH simula- 
tions used as input in SUNRISE. An order-of-magnitude increase 
in this resolution makes the face-on ultraviolet attenuation smaller 
and the edge-on near-infrared one larger (both by about 10 percent) 
for their model Sbc galaxy. J10 conclude that high-resolution sim- 
ulations are required to better explore the origin of these effects. 

2. 1.4 The S98 model iSilva et alll99&) 

The S98 model relies on the GRASIL radiative transfer code 
dSilva et aljfl998l ; iBressan et alj|2002l : IPanuzzo et ai]|2003h . This 
solves the radiative transfer equation in any location of model 
galaxies described by analyt ic distributions of dust an d stars, us- 
ing the iterative algorithm of iGranato & Danesel (fl994j). S98 com- 
pute in a consistent way the properties of stellar populations and 
the masses of interstellar gas and dust using an analytic model of 
galaxy evolution, which includes infall of pristine gas and the en- 
richment of the ISM by successive generations of stars. They as- 
sume that galaxies are composed of a stellar bulge and disc, pene- 
trated by a dust disc. For the bulge, they adopt a King profile of the 
form p(r) oc [l+(r/h,) 2 ]~ 312 truncated at log(r//?,) = 2.2, which dif- 
fers from the exponentially declining and de Vaucouleurs profiles 
adopted by P04 and T04, respectively. 

As P04 and T04, S98 adopt double exponentially declining 
light and density profiles to describe the stellar and dust discs, re- 
spectively (equation [6j, assuming that both discs have equal scale- 
lengths and scale-heights (Table[T](. We note that S98 also include 
a component of dust in stellar birth clouds, which attenuates the 
radiation from newly bora stars. As in the case of the T04 and 110 
models above, we ignore this component here and consider only 
the attenuation of starlight by dust in the diffuse ISM. 

We consider the attenuation curves ff M (ff) computed by S98 
for five model galaxies: an Sb and an Sc galaxy (the S-series) 
available from the GRASIL websit^l, and three models obtained 
from the spectral fits of late-type spir al galaxies (M 100 , M 5 1 and 
NGC 6946; the N-series) presented in lsilva et all dl998h . All mod- 
els were evolved for 13 Gyr. The Sb and Sc galaxies differ by their 
timescales for gas infall and efficiencies of star formation, which 
make the gas reservoir be depleted on a longer timescale in the Sc 
model than in the Sb one. The N-series galaxies are described as 



5 http: //adlibitum.oat.ts. astro. it/silva/grasil/modlib/ 
modlib.html. 
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bulge-less, with roughly similar timescales for gas infall and effi- 
ciencies of star formation. 



2.2 Relation between slope of the attenuation curve and dust 
attenuation optical depth 

We have described the different prescriptions of P04, T04, J 10 and 
S98 to calculate the attenuation optical depth ¥f M (ff) of dust in the 
ambient ISM at an angle 9 between the line of sight and the nor- 
mal to the equatorial plane of a galaxy. We now show that these 
prescriptions are consisted with one another and yield a quasi- 
universal relation between slope of the attenuation curve at any 
wavelength and V-band attenuation optical depth of the dust. To 
achieve this, we fit the results of the four models by expressing the 
dependence of f' ( SM (0) on A at fixed 9 as 



T> S >) = C<*) 



0.55 um 



(7) 



where the exponent of the power law is a linear function of wave- 
length, 



\9) = A(9) + B(9) 0*/um - 0.55) . 



(8) 



This parametrisation provides excellent fits to the radiative transfer 
calculations of P04, T04, J10 and S98. In practice, we select the 
best-fit A{9) and B(9) corresponding to every combination of input 
parameters of each of the P04, T04, J10 and S98 models. We focus 
below on the slope of the attenuation curve estimated in this way at 
A = 0.55 and 1.60 fim, which we refer to as and n™ M (9). 

Fig-Sshows n l y M (9) as a function of t [ y M (9) (top panels), and 
n^ M (9) as a function of t^ M (9) (bottom panels), as predicted by the 
P04, T04, J10 and S98 models over the full explored ranges in face- 
on dust optical depth and galaxy inclination. The left-hand panels 
show the attenuation of the stellar disc components in the P04, T04 
and S98 models along with the 7 SPH simulations (with insepara- 
ble dic+bulge components) of J10. The right-hand panels show the 
attenuation of the bulge component in the P04, T04 and S98 mod- 
els and the same composite J 10 models as in the left-hand panels. 
Different colours refer to different face-on dust optical depths. At 
fixed face-on dust optical depth, the same symbol plotted at dif- 
ferent abscissae shows the same model seen at different inclina- 
tions (see caption for detail). Also shown for reference as horizon- 
tal lines in all pan els in Fig.f4]are th e slopes o f the mean Milky-Way 
to'Donnelll994h. SMC and LMC |Pei l992h extinction curves and 
the lCalzettil 120011) attenuation curve. 

Fig. [4^ shows that, in all models, the slope riy (0) of the 
optical attenuation curve decreases (i.e. the curve becomes shal- 
lower) when the V-band attenuation optical depth f [ y M (9) increases. 
The attenuation curve is predicted to be steeper than the reference 
ICalzettil fcOOlh (Milky-Way) curve at f^ M (9) < 0.5 (f' v SM (e) < 
0.25) and shallower at larger f£ M (0). Differences in the disc ge- 
ometries adopted by P04, T04, J10 and S98 induce small vertical 
offsets between the different models (symbols), without affecting 
the general trend. It is worth stressing that the overlap of identical 
symbols with different colours in Fig.|4^ implies that galaxies with 
different dust content and inclination but analogous f l y M (9) have 
similar attenuation curves. To better illustrate this effect, we plot 
in Fig.[5]the attenuation curves of pure-disc models with different 
dust contents seen at different inclinations, for which the V-band 
attenuation optical depth fy M (0) differs by less than 0.05, as com- 
puted by P04, T04 and S98. We choose the widest possible dynamic 



range in 9 that can be sampled at fixed f'^ M (#) by combining pure- 
disc S98 models with different dust contents (i.e. from 48° for the 
Sb-disc model to 90° for the NGC 6946 model). Then, we select 
similar ranges in 9 from the model grids of P04 and T04. Fig. [5] 
confirms that models with different dust content seen at different 
inclination but with similar t l y U (6) have similar attenuation curves 
across the whole range from ultraviolet, to optical, to near-infrared 
wavelengths. This implies that the quasi-universal relation between 
n l y M (8) and f l y M {9) of Fig(4]is not purely the result of orientation 
effects. 

We interpret the relation between »I? (ff) and f If M (9) in Fig.|4^ 
as the combination of two effects, one dominating at low atten- 
uation optical depths, the other at large ones. At small attenua- 
tion optical depths, scattering is the likely reason for the steep- 
ness of the attenuation curve. Fig.[T]shows that, for the Milky- Way 
type dust in the P04, T04 and J10 models, the single-scattering 
albedo is roughly constant across the optical wavelength range 
3.6 < log(/i/A) < 3.9, while the asymmetry parameter of the scat- 
tering phase function decreases monotonically. Thus, scattering is 
more forward in the blue and more isotropic in the red. This implies 
that blue photons emitted along the equatorial plane of a galaxy will 
have less chances than red photons to be scattered away from the 
plane and escape from the galaxy before they are absorbed. To a 
face-on observer, this will appear as a steepening of the attenua- 
tion curve. We highlight that this effect is important only at small 
fy M (0) (corresponding also to small 9; see Fig.|4^), since at large 
M (0), photons are more likely to be absorbed even if they are first 
scattered, which reduces the imprint of scattering on the shape of 
the attenuation curve. Also, we note that the attenuation curve of 
the S98 model, in which scattering is assumed to be isotropic at 
all wavelengths (Fig. [TJ;), does not exhibit any steepening at low 
if M (6), reaching at most n l ^ M (9) « 1. We can illustrate the effect 
of non-isotropic scattering on the attenuation cu rve by appealing to 
the pu blicly available radiative transfer code of iMacLachlan et aU 
J201 1[), which relies on a Monte Carlo algorithm to solve the equa- 
tion of transfer for analyt ic distributions of s tars a nd dust. As in 
the P04 and T04 models, IMacLachlan et alj 1201 ll) adopt Milky- 
Way type dust and double exponentially declining profiles to de- 
scribe the distributions of stars and dust. Fig.|6]shows n l y M (9) ver- 
sus f y M (9) for two sets of disc calculations kindly made available to 
us by J. MacLachlan (private communication): a standard set with 
ranges in and 9 similar to those of the P04, T04 and J 10 models 
in Fig.|4^, in qualitative agreement with these models (open trian- 
gles in Fig.[6]l; and an identical set in which the single-scattering 
albedo has been set to zero (open circles). The results in Fig. [6] 
shows that scattering is the likely reason for the steepening of the 
attenuation curve at low r$ u (6) in Fig.|4^. 

At large attenuation optical depths, the shallowness of the at- 
tenuation curves in Fig. [4^ is a w ell-known signature of mixed 
distri butions of stars and dust (e.g. ICharlot & FalfeoOOt ICalzettil 
1 20011) : at any wavelength, the radiation reaching an observer from 
such distributions is dominated by stars located at an attenuation 
optical depth less than unity from the edge of the disc (which cor- 
responds to deeper physical locations for red photons than for blue 
photons). This causes the attenuation curve of the disc compo- 
nents) to be greyer than the input (Milky-Way) extinction curve 
at large f 1 y M (9) in the models of Fig. [4^. Fig. |4p shows that, in 
the P04, T04 and S98 models, the attenuation curve of the bulge 
component at small and large ry M (ff) behaves in a way similar to 
that of the disc component(s) in Fig. [4^, as expected from the fact 
that the same double exponentially declining dust disc(s) attenuate 
bulge and disc stars in these models (see Figs [2] and [3}- At fixed 
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Figure 4. (a) Slope of the attenuation curve «^ SM (0) at galaxy inclination 6 plotted against dust attenuation optical depth fy ((?) (equations [Vj£8j , at A ■■ 



0.55 



ujn, as predicted for the stellar disc components of the P04, T04 and S98 models and the (inseparable) disc+bulge components of the J10 SPH simulations. 
For the P04 and T04 models, different colours refer to different face-on optical depths r™^ 1 and tbx, respectively: 0.25 (cyan), 0.50 (blue), 1.0 (steel blue), 2.0 

' '" J v SM >e = 0.06), Gl (blue: (f™ 

0.48). For the S98 model, different colours 



(purple), 4.0 (orange) and 8.0 (red). For the J 10 model, different colours refer to different simulations: GO (cyan: {f^ M ), 
(steel blue: (tf M ) e = 0.14), G3 (purple: <f^ M ) 8 = 0.24), Sbc- (orange: (f I l/ SM > fl = 0.31), and Sbc+ (red: (tf M ) e 



refer to different simulated galaxies: Sb (cyan: (f y )g 



: 0.65), Sc (blue: (f™ 



■■ 0.97), M100 (steel blue: (ff M ) e - 

and NGC6946 (orange: (ff M ) 9 = 0.12). At fixed dust optical depth, the same symbol plotted at different (9) corresponds to the same model seen at 
different inclinations, 6 = °, 36 °, 50 °, 63 °, 74 °, 85 ° and 90 ° from left to right, (b) Same as (a), but for the bulge component of the P04, T04 and S98 
models [the J10 simulations are the same as in (a)], (c) Same as (a), but at A = 1.6 [im. (d) Same as (b), but at A = 1.6 [im [the J10 simulat ions are the same 
as in (c)]. For refer ence, in p anels (a) and (b), ho rizontal lin es indicate the slopes in the range 0.4 ^ /i/ u.m ^ 0-5 of the mean Milky-Way iO'Donnelj 1994, 
dotted line), SMC (Pei 1992, dashed line), LMC iPeill 19921 dot-dashed line) extinction curves and of the[Calzetti (2001, long-dashed line) attenuat ion curve. 
Horizontal lines in panels (c) and (d), indicate the slopes in the range 1 .4 ^ A/\im ^ 1.6 of the mean Milkv-Wav lPelll99l dotted line), SMC |Pelll99l 
dashed line), LMC (Pei 1992, dot-dashed line) extinction curves and of the Calzetti (2001, solid line) attenuation curve. The solid black line in (a) shows the 
mean relation given by equation )10K and the solid grey lines show the associated 25-percent dispersion [these lines are repeated in (b) for reference]. The 
solid black line in (c) shows the mean relation given by equations t9l- )l It , and the solid grey lines show the associated 25+10-percent dispersion [these lines 
are repeated in (d) for reference] . 



0.15), M51 (purple: <f I ™) 9 = 0.31) 



f'y M (#), the bulge attenuation curves predicted by the P04 and S98 
models in Fig. |4j5 are systematically shallower than that predicted 
by the T04 model. This is likely to be a consequence of the dif- 
ferent light profiles adopted for the stellar bulge in these models, 
which are much shallower in the P04 and S98 models than in the 
T04 model, as illustrated by Fig. [7] As we shall see in Section [3T41 
below, bulge stars in nearby star-forming galaxies account for typ- 
ically only a few percent of the total optical emission (Table 0. 
Thus, in practice, the global dust attenuation curve in these galax- 



ies should be well represented by that pertaining to the disc com- 
ponents) in Fig.|4ji. 

Figs[4j; and[4ji show the analog of Figs|4^ and|4p for the near- 
infrared attenuation curve. Again, all models in Figs |4j; appear to 
follow a quasi-universal relation between n^ M (d) and f^ M (0). The 
steepening and flattening of the attenuation curve at small and large 
f^ M (0), respectively, is qualitatively similar to the behaviour of the 



© 0000 RAS, MNRAS 000, 000-000 



8 J. Chevallard, S. Chariot, B. Wandelt, V. Wild 




\/fim 



Figure 5. Attenuation curves f' SM (0) of pure-disc models with different 
dust contents seen at different inclinations, for which the V-band attenua- 
tion optical depth t^ M (0) differ by less than 0.05. (a) For three P04 mod- 
els with tf u (g) * 0.55: (T™ d 



= (0.5,85°), (2.0,70°) and (8.0,35°). 
f M (0) a 0.45: (Tax.ff) = (0.5,85°), 
isMj-^-i ^ o 33' 



(b) For three T04 models with 

(2.0,75°) and (8.0,36°). (c) For three S98 models with r£ M (6) 
(disc model, 8) = (Sb disc, 48°), (M51,70°) and (NGC6946,90°). In each 
panel, the dashed, dotted and solid lines refer to the smallest, intermediate 
and largest inclination angles, respectively. 
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Figure 6. Slope of the attenuation curve n' SM (0) at galaxy inclination 8 
plotted against dust attenuation optical depth f 1 ™^) (equations |7J{8), at 
A = 0.55 um, as predicted by the MacLachlan et al. ( 20 11) radiative transfer 
code for disc galaxies. Open triangles show the results for standard calcu- 
lations, while open circles show identical calculations in which the single- 
scattering albedo has been set to zero (MacLachlan, private communica- 
tion). In each case, the different points correspond to galaxies with different 
face-on V-band optical depths, t™'^ = 0.1, 0.5 and 2.0, and inclinations, 
from 8 = ° to 90 ° in steps of 10 ° . 



optical attenuation curve at small and large f [ y M (S). In Fig. |4jl, 
the behaviour of the bulge attenuation curve at near-infrared wave- 
lengths is also qualitatively similar to that found in Fig.[4p at optical 
wavelengths. Again, the different trends followed by the P04 and 
S98 models on the one hand, and the T04 model on the other hand, 
are likely to result from the different light profiles adopted for bulge 
stars in these models (Fig. [7}. As before, we note that this discrep- 
ancy should have a negligible implication for measurements of the 



o 



•8 0.5 




Figure 7. Light profiles adopted for the stellar bulge component in the P04 
[exponentially declining profile: p(r) oc exp(-r//u), where h„ is the bulge 
scale length], T04 jde Vaucouleurs profile: p(r) oc exp[(r//i,)'' 4 ]} and S98 
(King profile: p(r) oc [1 + (r/A«) 2 ] _3 ' 2 | models. 



global attenuation curve in nearby star-forming galaxies, which is 
expected to be dominated by the disc component(s) (Section|33J. 

For practical applications, it is of interest to fit the quasi- 
universal relation between slope of the attenuation curve and at- 
tenuation optical of the dust in the ambient ISM in Figs|4^ and|4j; 
by means of a simple analytic function. Adopting for this purpose 



the simplified notation r v 
tion l[8} as 



(0), we find that rewriting equa- 



nf M (T V ) = nlfVv) + b(r v ) Ol/nm - 0.55) , 



(9) 



with 



and 



2.8 



1 +3Vt> 



(±25 percent) (10) 



b(T V ) = 0.3 - 0.05 t v 



(±10 percent) , 



(11) 



provides an adequate representation of the relations both between 
4 SM ((9) and f® M (0) and between n^ M ((9) and f™((9), as shown by 
the black solid lines in Figs|4ji and|4j; (grey solid lines indicate the 
dispersion around these mean relations). We emphasise that equa- 
tions (7) and d9t — d 1 It summarise the results of four sophisticated 
radiative transfer models including realistic spatial distributions 
of stars and dust. These equations can be easily implemented in 
any spectral analysis of structurally unresolved star-forming galax- 
ies to include in a simple yet physically consistent way for the 
dependen ce of dust attenuat ion on inclination (as exemplified re- 
cently bv lPacifici et"aT]|2012l) . The fact that the P04, T04, J10 and 
S98 models all rely on standard optical properties of Milky Way- 
type dust should be kept in mind, but we do not expect this to 
strongly limit applications of equation (|9} to spectral studies of ex- 
ternal galaxies. In fact, the relation in Fig.|4^ extends over a much 
wider range in n\ f M (0) than the d ispersion between the Milky-Way, 
LMC, SMC and lCalzettil d200ll) curves. Also, the optical proper- 
ties of Milky Way- and SMC-type dust are similar at optical wave- 
lengths, exhibiting larger differ ences in the ul traviolet (see fig. 10 of 
iGordon et al.l2003l and fig. 4 o flDrainel2003l) . This suggests that, in 
normal star-forming (disc) galaxies, changes in the shape of the op- 
tical attenuation curve induced by geometry and orientation effects 
are likely to dominate over those induced by differences in the op- 
tical properties of dust grains. Finally, we note that the dispersion 
about the mean relations in equations dlOt and dl It accounts to 
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some extent for potential systematic uncertainties in these relations 
arising from the different assumptions of the P04, T04, J10 and 
S98 models (i.e. different dust properties, different scale-lengths 
and scale-heights of the stellar and dust discs; see Table[TJ. A more 
thorough investigation of the assumptions behind radiative-transfer 
models and their implications for the systematic uncertainties af- 
fecting the results presented in this paper would require running 
the original models for different values of their main adjustable pa- 
rameters. This is beyond the scope of the present study. 

2.3 Relation between dust attenuation optical depth and 
galaxy inclination 

In the previous section, we found that the sophisticated radiative 
transfer models of P04, T04, J 10 and S98 all predict a similar de- 
pendence of the slope of the attenuation curve on the dust atten- 
uation optical depth of the diffuse ISM in star-forming galaxies. 
Here, we show that, at fixed dust content, these models also pre- 
dict a similar dependence of t l f M (8) on galaxy inclination 6. We 
demonstrate this by using the more versatile T04 model, which in- 
cludes a thin-disc, a thick-disc and a bulge components, to repro- 
duce the attenuation curves of different types of galaxies computed 
with the P04 and J10 models (we do not consider here the S98 
model, in which the approximation of isotropic scattering makes 
the attenuation curves much shallower than in the other models at 
low optical depths, i.e. at face-on inclinations, which would bias 
the comparisons). To proceed, we explore a wide collection of at- 
tenuation curves generated using the T04 model by considering full 
ranges of central face-on 6-band optical depth t b± (Section [2.1.2t 
and of relative intensities fj of the thin-disc, thick-disc and bulge 
components. In practice, for a given star formation law i/f(f), we 
compute £j from equations ((3) and ((4) by relating stars in dif- 
ferent age ranges to different geometric components of the T04 
model. This is achieved by introducing three age parameters, t thin , 
Wk and tbuige, such that tthin ^ tmct ^ Wge- We associate stars 
younger than t thin with the thin-disc component, those in the age 
range t lhin < t sC t thick with the thick-disc component and those 
in the age range tthick < t ^ tbuige with the bulge component (we 
consider the ranges in t b± , t thin , t thick and t bll i ec listed in TableO. To 
compute the spectral evolution S ,\ [f ,Z(t-t')] of stellar popu lations 
in equation ©, we appeal to the lBruzual & Chariot! J2003) stellar 
population synthesis code. 

We focus on the dependence of the attenuation optical 
depth f^ M (ff) on inclination at the effective wavelengths of the 
ugrizYJH photometric bands (A cff = 3546, 4670, 6156, 7471, 8918, 
10305, 12483, and 16313 A, respectively). To efficiently explore 
the T04 parameter space and find the combination of t b ±, tthin, tthick 
and t bu i ge that can best reproduce calculations of ff u (6) by P04 and 
J10, we appeal to a Markov Chain Monte Carlo (MCMC) algorithm 
of the type described in Section [331 below (to which we refer for 
detail; see also Appendix [Cl. In Fig. [8] we show examples of two 
combinations of T04 model parameters selected in this way to re- 
produce P04 predictions for two pure-disc models with face-on V- 
band dust optical depth t™^' = 1 and 4, respectively (Section |2.1.U . 



mid _ i 
l V± ~ 1 



In practice, we fit the P04 pure-disc model with the T04 thin- and 
thick-disc components, excluding the T04 bulge component (i.e. 
we fix tthick = tbuige)' F° r simplicity, since the P04 model does not 
include any star formation history, we have adopted a constant star 
formation law, if/(t) = const, to compute the relative intensities £j 
from equation ((4} in the T04 model. The excellent agreement be- 
tween the stars and circles in all panels of Fig.[8]demonstrates that 
two different combinations of the T04 thin- and thick-disc compo- 
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Figure 8. Dust attenuation optical depth fJ 1 SM (0) at galaxy inclination 8, in 
units of the angle-averaged value (f' SM )g, plotted against 1 - cos 8, in the 



ugrizYJH photometric bands, as indicated. Top: reproduction of a pure-disc 
P04 model with r^ d = 1 (stars) by a T04 model with t b ± = 2.2 and age 
parameters t t hin = 0.01 Gyr and t t hick = tbuige = 6.2 Gyr (open circles). 
Bottom: reproduction of a pure-disc P04 model with rjjjj = 4 (stars) by 
a T04 model with t b ± = 5.7 and age parameters t t hi n = 0.003 Gyr and 
tthick = tbuige = 9.7 Gyr (open circles). See Section |2~3l for detail. 



nents can reproduce the dependence of f' SM (0) on 6 predicted at all 
wavelengths by P04 bulge-less models with = 1 and 4. We 
note that the T04 best-fit parameters (listed in the figure caption) 
depend on the adopted star formation law if/(t), which cannot be 
constrained independently using the P04 attenuation curves. 

In Fig. [9] we present fits of the Sbc+ and G3 models of J10 
obtained in a similar way by optimising the combination of T04 
model components. We choose these galaxies because they are 
the largest ones simulated by J10 in their respective mass ranges 
(Section [2.1.3l >. which should minimise the influence of numerical 
resolution on the results. To compute the relative intensities £j of 
the different geometric components of the T04 model in this case 
(equation^, we have adopted the original star formation histories 
of the J10 simulations. Fig.[9]shows that the fit of the Sbc+ model is 
overall excellent, despite some slight deviations in the near-infrared 
Y, J and H bands. The deviations in these bands are more pro- 
nounced for the fit of the G3 model. In fact, the G3 model reaches 
large negative attenuation optical depths at low inclination in the 
Y, J and H bands. This is again a signature of the isotropic scat- 
tering of infrared photons emitted along the equatorial plane of a 
galaxy, which can be scattered away from the plane and boost the 
emission in the face-on direction. The effect is more pronounced 
in the G3 model than in the Sbc+ model, possibly because of the 
lower dust attenuation optical depth of the former «f I l f M ) e = 0.24) 

0.48; see the discussion in Sec- 



relative to the latter {(f^ M ) e 
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Figure 9. Dust attenuation optical depth f' SM (0) at galaxy inclination 8, 
in units of the angle-averaged value (f' ( SM )e, plotted against 1 - cos0, in 
the ugrizYJH photometric bands, as indicated. Top: reproduction of the 
Sbc+ model of J10 (stars) by a T04 model with tb± = 3.0 and age pa- 
rameters t t hin = 4.7 Gyr, tthick = 10.8 Gyr and t^ige = 11.3Gyr (open cir- 
cles). Bottom: reproduction of the G3 model of J10 (stars) by a T04 model 
with Tgx = 1.8 and age parameters t t hi n = 0.5 Gyr, t t hick = 11.3 Gyr and 
tbuige = 11.5 Gyr (open circles). See Section l231 for detail. 



tion l2.2l above). We note that the T04 model also reaches negative 
infrared attenuation optical depths at = 0, but only of the order of 
a few hundredths of magnitude. This discrepancy between the in- 
frared predictions of the T04 and J 10 models has a negligible effect 
on galaxy spectral analyses, since the infrared attenuation optical 
depths are generally very small (for reference, the angle-averaged 
/f-band attenuation optical depths of the Sbc+ and G3 models are 
{r l H M (d))e = 0.11 and 0.03, respectively). 

In summary, we find that the P04, T04 and J 10 models all pre- 
dict a roughly universal relation between attenuation optical depth 
of the dust and slope of the optical attenuation curve in the diffuse 
ISM of star-forming galaxies. We also find that, at fixed dust con- 
tent, the three models predict a similar dependence of fJ, SM (0) on 
galaxy inclination 6 in the ugrizYJH photometric bands. We could 
demonstrate this because the versatility of the T04 model allows 
us to reproduce calculations of ff M (6) by P04 and J10 at all these 
wavelengths for different types of galaxies. 



as a function of both inclination and star formation history in a sam- 
ple of nearby star-forming galaxies. We start with a brief descrip- 
tion of the observational sample. Then, we build a comprehensive 
library of model spectral energy distributions encompassing wide 
ranges of stellar and dust parameters to interpret the Ha/H/3 ratio 
and ugrizYJH photometry of the observed galaxies. This allows us 
to derive original constraints on the content and spatial distribution 
of dust from the integrated spectral properties of these galaxies. 



3.1 Observational sample 

3.1.1 Sample selection and basic parameters 

We conside r the large s ample of nearby star-forming galaxies as- 
sembled bv lWild etallfcoilbl. hereafter Wll) by combining op- 
tical spectros copy and ugriz photom etry from the SDSS Data Re- 
lease 7 (DR7, Abazai ian et al J2009I) with near-infrared YJHK pho- 
tometry from the seventh data release of the UKIR T Infrared Deep 
Sky Survey-Large Area Survey (UKIDSS-LAS, lLawrence et al.l 
120071) . To build this sample, Wll first select SDSS galaxies with 
high-quality spectra (S/N > 6 per pixel in the g band) within the 
formal survey limits (14.5 < m, < XI 11, where m, is the r- 
band Petrosiarf] magnitude corrected for Galactic extinction) and 
with physical-parameter entries in the SDSS-MPA value-added cat- 
alogQ Potential host s of act ive galactic nuc lei are removed usini 
thelKauffmann et al. ( 2003c) criterion in the lBaldwin et al.l dl98ll 
diagram defined by the [N ii]/Hcj' and [O m]/HyS emission-line ratios 
(with the requirement that all four lines be measured with S/N > 3, 
after rescaling of the noise as recommended on the SDSS-MPA 
website). Wll further exclude g alaxies with pos t-starburst spectral 
features using the method of Wild et al.l d2007t) . since the rapidly 
changing spectral shapes of these galaxies can affect measurement 
of the attenuation curve through the pair-matching technique (see 
Section l3T2l below). Wl 1 cross-match the SDSS and UKIDSS cat- 
alogs by locating all objects within 2 arcsec of each other, selecting 
the nearest neighbour in case of multiple matches. They exclude 
galaxies for which UKIDSS photometr y has been de- blended, as 
these suffer from flux calibration issues jffill et alj20icb . The final 
sample consists of 22,902 galaxies at a median redshift of 0.07. 

The median r-band Petrosian radius of these galaxies is 
5.1 arcsec, i.e. larger than the 3 arcsec-diameter aperture of the 
SDSS spectroscopic fibres. Since we are interested in physical 
quantities derived from both spectroscopy and photometry, we do 
not use Petrosian magnitudes here but adopt instead fibre mag- 
nitudes measured in a 3 arcsec-diameter aperture from the SDSS 
catalog. We combine this with infrared photometry measured in a 
matching aperture of 2.83 arcsec diameter from the UKIDSS cata- 
log (see Wl 1 for more detail). 

Basic physical parameters are available for these galax- 
ies from published SDSS analyses. This includes both total 
(M») and fibre (M^ b ) stellar masses, measured from Bayesian 
fits to the five-band SDSS photometry using a comprehen- 
sive library o f mod el galaxy spectra similar to that used in 
iGallazzi et all d2005L accounting for dust attenuation as pre- 



3 INSIGHTS INTO THE CONTENT AND SPATIAL 

DISTRIBUTION OF DUST IN OBSERVED GALAXIES 

In this section, we exploit the general properties of dust attenuation 
models outlined in the previous section to investigate the physical 
origin of observed systematic changes of the dust attenuation curve 



6 This corresponds to the total flux within a circular aperture of twice the 
Petrosian radius, where the Petrosian radius is the largest radius at which 
the local r-band surface brightness is at least one-fifth the mean surface 
brightness interior to that radius. 

7 http: //www.mpa-garching .mpg .de/SDSS/DR7/raw_data.html 
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Figure 10. Distributions of physical properties of the low-//* (solid his- 
tograms) and high-//, (dashed histograms) nearby star-forming galaxies in 
the Wll sample described in Section [3. 1.1 1 (a) Gas-phase oxygen abun- 
dance, 12 + log(0/H). (b) Specific star formation rate, i//s- (c) Stellar mass, 
logM». (d) Radius containing 90 percent of the r-band Petrosian flux, Rgg. 
The distributions of 12 + log(0/H) and tf/s are used as priors in the library 
of model spectral energy distributions described in Section POl 



scribed by Chariot & Fall 2OOO)0and gas-phase oxygen abundance 
[12 + log(0/H); iTremonti et alj|2004[ correct ed for dust attenua- 
tion using the model of ICharlot & Falllboodl and star formation 
rate (if/; we adopt here the fibre estim ates of Wl 1, correc ted for dust 
attenuation using the prescription of Iwild et alj|201 lah . measured 
from fibre spectroscopy. From these quantities, we compute the 
specific star formation rate within the fibre aperture, i/fs = (fr/Mj b , 
and the stellar surface mass density, fi, = M t /(2nR^ _), where 
i?5o,c is the radius encompassing 50 percent of the z-band Pet- 
rosian flux. Wll differentiate galaxies according to p.,. Following 
Kauff mann et "all( l2003j) . they distinguish bulge-less galaxies, with 
p, < 3 x 10 s Mqpc~ 2 (7451 galaxies), from galaxies with a bulge, 
with p., > 3 x 10 s M pc~ 2 (15 451 galaxies). We show in Fig.fTOl 
the distributions of 12 + log(0/H), t/r$ and M, for these two sub- 
samples, along with the distribution of the radius enclosing 90 per- 
cent of the total Petrosian r-band flux, R 90 . The low-//, galaxies ex- 
tend to significantly lower gas-phase oxygen abundances and stel- 
lar masses than their high-/j, counterparts, although in both cases, 
similar ranges in specific star formation rate and size are probed. 
Also available from the SDSS catalog for these galaxies is the ra- 
tio of minor- to major-axis length (bid) measured from exponential 
fits to the SDSS r-band image. We convert this ratio into inclina- 
tion angle u sing the standard formula for an oblate spheroid (e.g. 
iGuthriell 19921) 



cos 6 -- 



(bla? 



. 2 

1 ~% 



(12) 



where q is the intrinsic flattening of the ellipsoid representing the 
galaxy (i.e. the axis ratio when seen edge-on). We adopt q = 0.14 
and 0.20 fo r the low- and high-//, samples, respectively llRvdenl 
l2004l2006h . 



http: //www.mpa-garching .mpg .de/SDSS/DR7/Data/ 
stellarmass .html 



3. 1.2 Analysis o AWildet al\ \2011b\) 

Wll use the above sample to explore the dependence of the shape 
of the optical and near-infrared attenuation curves on galaxy incli- 
nation and specific star formation rate. They achieve this by com- 
paring the observed spectral energy distributions of pairs of galax- 
ies with presumably similar star formation and chemical enrich- 
ment histories, but different dust content (as measured from the 
observed Ha/H/3 ratio), at fixed i nclin ation. This is similar to the 
approach introduced by ICalzettil d200ll) to constrain the attenua- 
tion curve in nearby starburst galaxies, only generalised to galaxies 
with more complex star formation histories. In essence, Wll first 
identify pairs of galaxies with comparable gas-phase metallicity, 
specific star formation rate, axis ratio (inclination) and redshift but 
different dust content (replacing the match in gas-phase metallicity 
by one in stellar mass or surface mass density does not alter the con- 
clusions). They normalise each spectral energy distribution to unit 
flux in the K band (A?g = 22010 A) to remove the dependence on 
galaxy mass. Then, the flux ratios of the two galaxies in the SDSS 
ugriz and UKIDSS YJH bandpasses provide an estimate of the dust 
attenuation curve. Wll compute average flux ratios in bins of fixed 
inclination b/a and specific star formation rate i/r s to measure the 
dependence of the optical and near-infrared attenuation curves on 
these parameters. 

The conclusions of Section [2~2l above allow us to reconsider 
the results of Wl 1 from a new perspective. Mean atte nuation curves 
derive d from differential photometry, such as those of Calze tti et al.l 
d 19941) and Wll, reflect the true attenuation curve of a class of 
galaxies only if the shape of this curve does not depend on the 
Ha/Hfl ratio. However, as shown by Fig. [4] current sophisticated 
models of radiative transfer predict that the shape of the dust attenu- 
ation curve should change with dust content in galaxies with similar 
relative distributions of stars and dust. Thus, differential photome- 
try of pair-matched galaxies with similar intrinsic spectral energy 
distributions and different Ha/H/3 ratio should trace not only the 
dust attenuation curve, but also changes in the shape of this curve 
with absolute dust content. The framework developed in Section [2] 
allows us to account for these two contributions and better charac- 
terise the dependence of the attenuation curve on galaxy geometry 
and inclination. We note that Wll also analyse the dependence of 
the attenuation curve on specific star formation rate. We do not in- 
vestigate this dependence here because the radiative transfer mod- 
els described in Section |2T| do not include any physical relation 
between attenuation by dust and star formation rate. 



3.1.3 Correction of systematic biases in the dependence of dust 
attenuation signatures on inclination 

The first step toward characterising the dependence of the atten- 
uation curve on galaxy geometry and inclination is to accurately 
measure the changes of observed spectral properties with axis ra- 
tio. We consider here the broadband spectral energy distribution, 
which traces the dependence of attenuation on wavelength, and 
the Ho-/H/3 ratio, which quantifies the absolute attenuation opti- 
cal depth. The Ha-/H/J ratio reflects the attenuation of young stars 
by dust both in stellar birth clouds and in the ambient ISM (this 
ratio does not reveal the potential minor fraction of dust in the 
ionised gas, which can abso rb photons before they ionise hydro- 
gen; e.g. Icharlot&FaUE oOO). We measure the attenuation relative 
to the dust-free case B recombination ratio Ha/H/3 ~ 2.86 (ap- 
propriate for electron densities n e < 10 4 ctrr 3 and temperatures 
T e as 10 000 K), the difference between the dust attenuation optical 
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depths at H/3 and Ho- being given by 



THf3(0) ~ T Ha (d) = In 



Hg/HM 
2.86 j 



(13) 



To proceed, we first divide the sample of Section 13.1.1 1 in 10 bins 
of equal axis ratio b/a. By analogy with Wll, we normalise the 
ugrizY J HK spectral energy distribution of each galaxy to unit flux 
in the K band. Then, we compute the mean spectral energy distri- 
bution and mean Ha /H/3 ratio in each bin of b/a. By comparing the 
mean spectral properties of galaxies in any bin to those in the face- 
one bin (b/a as 1; corresponding to the smallest attenuation), we 
can quantify the dependence of attenuation on galaxy orientation. 

In Figs [TT] and [12] we show the Ha /H/3 ratio and differen- 
tial ugrizY JH magnitudes, Am = m(b/a) — m(face-on), obtained 
in this way as a function of b/a, for the low- and high-^/, samples, 
respectively (open squares). As expected, when b/a decreases, i.e. 
as galaxies are seen more edge-on, the Ha/Hfl ratio rises, the ul- 
traviolet and optical flux drop significantly (Am > 0) while the 
near-infrared flux remains roughly constant (Am ss 0). To discrim- 
inate the true influence of inclination on these changes, we must 
ascertain that other physical parameters, which correlate with dust 
content, do not also systematically vary with b/a. For this pur- 
pose, we show in Figs 113b and 114b the difference between bin- 
averaged and sample-averaged values of several physical parame- 
ters as a function of b/a, for the low- and high-yu, samples, respec- 
tively. The parameters we consider are: gas-phase oxygen abun- 
dance, 12 + log(0/H); specific star formation rate, i/rs; redshift, z; 
and radius enclosing 90 percent of the total Petrosian r-band flux, 
Rgo- 

Figs|13b and!14b reveal important biases in the physical prop- 
erties of galaxies seen at different inclinations. In both the low- and 
high-/i, samples, galaxies seen more edge-on tend to have on av- 
erage systematically lower 12 + log(0/H), if/ s and z and larger R go 
than galaxies seen more edge-on. This is likely a consequence of 
the high signal-to-noise ratio required to detect the [N n], Ha, [O in] 
and Hf} emission lines in Section [3. 1 . 1 | (al though part of the trends 
could also result from measurement biases). In particular, since the 
dust attenuation optical depth increases with inclination, and since 
the most metal- rich galaxies also contain the most dust (see fig. 6 of 
iTremonti et"aT]|2004l and fig. 6 of iBrinchmann et al.ll2004T) . at high 
inclination, nebular emission lines of metal-rich galaxies are more 
difficult to detect than those of metal-poor galaxies. The fact that 
the specific star formation r ate correlates with gas - and hence dust- 
mass fraction (see Fig. 6 of Ida Cunha et al . 2010) produces a sim- 
ilar bias against highly-inclined, actively star-forming galaxies in 
Figsll3b and[l4] A bias against small Rgo at high inclination arises 
from a positive correlation between Rg and surface brightness for 
the galaxies in our sample, and possibly from the blurring of the 
shapes of the smallest galaxies by the observational point spread 
function. This blurring may also be responsible for the bias against 
highly inclined galaxies at large z, together with cosmological dim- 
ming. These different biases are substantial, reaching roughly 40 
(20) percent of the sample-averaged value in 12 + log(0/H), 30 
(60) percent in i//$, 40 (40) percent in z and 50 (70) in Rg for low- 
/j* (high-p.) galaxies. Such biases must be corrected to properly 
investigate the dependence of dust attenuation on galaxy inclina- 
tion. 

In the absence of a reliable model to describe complex corre- 
lations between 12 + log(0/H), i//$, z and Rgo, we can reduce the 
strength of the biases in Figs I13fa and [T4b using the data alone. 
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Figure 11. Top: Mean Ha/H/3 ratio in bins of different galaxy inclination, 
from face-on (1 - b/a a 0) to edge-on (1 - b/a al), for the low-/i, galaxies 
in the Wll sample described in Section 13.1.1 1 Square correspond to di- 
rect observations, while circles with error bars show the effect of correcting 
these observations for systematic biases using the method of importance 
sampling described in Section |3.1.3l and Appendix 151 Bottom: difference 
between the mean magnitude of galaxies in bins of different inclination and 
that of face-on galaxies, at fixed K band flux (i.e. fixed stellar mass), in the 
ugrizY JH bandpasses, as indicated. The squares and circles have the same 
meaning as in Top. 
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Figure 12. Same as Fig. II II but for the high-/j, galaxies in the Wll sample 
described in Section |3. 1.11 
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Figure 13. Relative difference between the mean physical property of 
galaxies in bins of different inclination, from face-on (1 - b/a k 0) to edge- 
on (1 - b/a » 1), and those of face-on galaxies, for the low-p* galaxies 
in the Wll sample described in Section |3 . 1 . 1 U dark green: 12 + log(0/H); 
blue: i//$; red: z; dark grey: Rgo)- (a) As measured in the original Wl 1 sam- 
ple, (b) After correcting these measurements for systematic biases using the 
method of importance sampling described in Section l3.1.3l and AppendixlBl 
For reference, the dashed horizontal lines show the effect of a ±10 percent 
bias. 



Changes in the mean value of a given parameter from a bin to an- 
other in these figures reflect differences in the underlying distri- 
butions of that parameter at different b/a. The combined changes 
in mean 12 + log(0/H), z and R 90 therefore reflect differences 
in the joint distributions of these parameters. To remove the in- 
fluence of these differences on the relation between dust atten- 
uation and galaxy inclination, we require a joint distribution of 
[12 + log(0/H), z, R90] common to all b/a bins in Figsll3b and 
114b . This can be extracted from the original sample by appealing to 
the formalism of importance sampling, a well known Monte Carlo 
tech nique used to solve, for instance, multi-dimensional integrals 
fe.g. lRobert & Casellal2004h . This method requires two conditions: 
that the joint probability distribution of [12 + log(0/H), tf/ s , z, R] 
be evaluated in each bin (this implies computing probability dis- 
tributions from a finite sampling of observables); and that the in- 
tersection of the joint probability distributions of different bins be 
non-zero. In Appendix [B] we describe in detail our implementa- 
tion of this technique to derive a joint probability distribution of 
[12 + log(0/H), i/rs, z, R90] common to all bins of b/a. The mean 
spectral properties of galaxies obeying this distribution are shown 
as a function of galaxy inclination by the open circles in Figs [TT1 
and!12l The associated error bars reflect the standard error of the 
mean. Figsll3b and 1 14b show that, after correction, any systematic 
change of a parameter with inclination amounts to at most about 10 
percent of the face-on value of that parameter, for both the low- and 
high-p, galaxy samples. 

Figs [TJJ and [12] show that the corrections of systematic biases 
in the physical properties of galaxies as a function of inclination 
have a larger impact on the Ha/Hj8 ratio and differential ugrizYJH 
magnitudes of low-yu, galaxies than on those of their high-/i„ coun- 
terparts, especially at high inclin ation. To investig ate the origin of 
this result, we plot in Fig. [B]the|Sfi earm artl ll 19041) rank correlation 
coefficient for the relations between Ho-/H/3 ratio and the four pa- 
rameters 12 + log(0/H), i/'s, z and i? 9 o, in each bin of b/a, for low- 
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Figure 14. Same as Fig. 1131 but for the high-/i» galaxies in the Wll sample 
described in Section |3. 1.11 
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Figure IS. Spearmar] J 19041) rank correlation coefficient for the relations 
between Ha/H/3 ratio and different physical properties of low-,u» (dia- 
monds) and high-jU* (triangles) galaxies in the Wll sample described in 
Section [3. 1.1 1 in bins of different inclination, from face-on (1 - b/a 1 0) 
to edge-on (1 - b/a as 1). (a) 12 + log(0/H); (b) )fr s ; (c) z and (d) R 90 . 
In each panel, filled and open symbols indicate a significance level of the 
correlation larger and smaller than 99.99 percent, respectively. 



and high-yu, galaxies. Figs 1 15b and c show that, for low-yu, galax- 
ies (diamonds), 12 + log(0/H) and z are strongly correlated with 
HQ'/Hyt> ratio at all inclinations. These are also the quantities ex- 
hibiting the largest biases at all inclinations in Fig. 113b. For high-^/, 
galaxies (triangles), Figs|15b-c show that 12 + log(0/H) and z are 
less strongly correlated with Ha/HfS ratio than for low-yU„ galaxies, 
while i^s is much more strongly correlated and R90 equally weakly 
correlated but in the opposite sense. The complex interplay between 
these different trends is likely to be the reason for the different cor- 
rections of the relation between Ha/HfS ratio and inclination per- 
taining to low- and high-/j, galaxies in Figs [TT] and [12] We can- 
not produce the equivalent of Fig.Q3]for the differential ugrizYJH 
magnitudes, because in any given bin of b/a, we compute only one 
measure of Am to average out the intrinsic variations in the spectral 
energy distributions of all galaxies in that bin (hence no correlation 
coefficient can be computed within a bin). 
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3.2 Library of model spectral energy distributions 

The dust attenuation models presented in Section[2]provide an ideal 
framework to interpret the observed dependence of the Ha/HyS ra- 
tio and ugrizYJH attenuation curve on galaxy geometry and in- 
clination (Figs [TT] and 1121 ). The four types of models we investi- 
gated in Section [2~T1 all predict similar changes of the shape of the 
attenuation curve in the ambient ISM as a function of V-band at- 
tenuation optical depth f l y U (9) (Section |2,2> . and of the attenua- 
tion optical depth f l ^ M (9) in the ugrizYJH bands as a function of 
galaxy inclination 6 (Section [2.3l l. To explore the physical origin 
of the observed trends in Figs II II and 1121 we appeal to the most 
versatile of these models, by T04, in which galaxies are described 
as the sum of three stellar components (a thin and a thick disc, 
and a bulge) and two dust components (associated with the stellar 
discs; see Section [2. 1.2b . In the next paragraphs, we combine the 
T04 model with a hierarchical galaxy formation model, a stellar 
population synthesis code and a prescription for dust attenuation in 
stellar birth clouds to build a library of spectral energy distributions 
tailored to the analysis of the sample of Section [3~Tl 



formly in redshift in the interval 0-1 (0-0.5) and rejecting galax- 
ies that formed less than 15 (10) percent of their stars in the last 
2.5 Gyr (parameter / SFH of lPacifici et alj20~12l) provides model dis- 
tributions in the D n 4000 - HSa diagram in qualitative agreement 
with the observations of low-yu, (high-yu.) galaxies. This is illus- 
trated by Figs|16fa and 1 16b for the low-yu„ and high-//, galaxy sam- 
ples, respec tively. 

As in IPacifici et al] d2012l) . we combine the two libraries 
of star formation and chemical en richment histories generat ed in 
this way with the latest version of iBruzual & Chariot! Hoolh stel- 
lar population synthesis code (Chariot & Bruzual, in preparation) 
to compute the spectral energy distribution of every galaxy at 
wavelengths between 91 A and 160/um . This code incorporates 
a new library of observed ste llar spectra at optical wavelengths 
(Sanchez-Blazquez et a i] |2006t) and new prescriptions for the evo- 
lution of stars less massive than 20 M© dBerteili et all 2008. 2009) 
and of thermally pulsing asymptotic giant branch (TP-AGB) stars 
(Marigo et al. 20 081). We adopt h ere the Galactic-disc stellar initial 
mass functionof lchabrieij d2003l) . 



3.2.1 Spectral energy distributions of stellar populations 

We fo llow the approach recently developed by IPacifici et alj 
d2012l) and appeal to state-of-the-art models of star-formation and 
chemical-enrichment histories to build a broad library of as realistic 
as possible spectral energy distributions of galaxies. We start from 
a catalog of 4000 star formation and chemical enrichment his- 
tories produced by IPacifici et alj d2012l) through the semi-analy tic 
post-treatment (based on recipes by foe Lucia & Blaizotll2007f) o f 
the Millennium cosmological simul ation of Springel et alj d2005l) . 
As discussed by Paci fici et alj d2012l see their section 2.1), star for- 
mation histories derived in this way are appropriate mainly to de- 
scribe the properties of nearby SDSS galaxies more massive than 
5 x 10 9 /t'Mq. The parameter space must be broadened to interpret 
observations of galaxies drawn from other samples and in other 
mass ranges, in particular by resampling the distributions of galaxy 
evolutionary stages (i.e. the age at which a galaxy is looked at) 
and current galaxy properties (i.e. the mean star formation rate and 
metallicity of new stars over a period of 10 Myr before the galaxy 
is looked at). Such resampling is compatible with the stochastic 
nature of star formation and ch emical evolution in hierarchical sce- 
narios of galaxy formation (see lPacifici et al.ll2012l for detail). 

From the 40 000 star formation and chemical enrichment his- 
tories mentioned above, we therefore produce two smaller catalogs 
of 5000 galaxies each, appropriate to analyse the observed low- 
and high-yu„ samples of Section [3TTl (we have checked that generat- 
ing twice-larger catalogs does not affect our conclusions). In prac- 
tice, we extract models one by one from the original library and 
randomly redraw evolutionary stage, specific star formation rate 
and metallicity of newly-bom stars from prior distributions appro- 
priate for the low- and high-yu, samples. Specifically, we take the 
prior distributions of specific star formation rate and metallicity of 
newly-born stars to be the distributions of if/ s and 12 + log(0/H) 
of the observed galaxies in both samples (Fig.[l0]of Section l3.1.1t . 
To redraw evolutionary stage, we appeal to the standard starburst- 
diagnostic diagram defined by the strengths of the 4000 A break 
and th e H<5a stellar absorption-line indices (e.g. iKauffmann et alj 
l2003al) . This diagram reflects the balance between newly bora stars 
with low D n 4000 and US a, intermediate-age stars with large US a 
and older stars with large D n 4000. After some experimentation, we 
find that resampling the age at which the galaxy is looked at uni- 



3.2.2 Inclusion of dust modelling 

We adopt the versatile radiative transfer code of T04 to describe 
the attenuation of starlight by dust in the ambient (i.e. diffuse) 
ISM (Section |2.1.2t . We must also account for the enhanced at- 
tenuation of newly born stars in their parent molecular clouds, 
whi ch dissipate on a timescale typic ally of the order of 10 Myr 
(e.g. lMurrav et al . 2010; Murrav 201 1 ). We achieve this by combin- 
ing the T04 m odel with the simple tw o-component, angle-averaged 
dust model of ICharlot & Falll d2000t) . which accounts for the dif- 
ferent attenuation affecting young and old stars in galaxies. We 
introduce a dependen ce on viewing angle in the prescription of 
Chariot & Fall ( 2000) and write the total attenuation optical depth 
at galaxy inclination 9 affecting the radiation from stars of age / as 



T7W) 



T ',1 

\e,t) 



(9,1) 



for f< 10 Myr, 
for t > 10 Myr. 



(14) 



Here f^ c is the attenuation optical depth, assumed isotropic, of the 
stellar birth clouds, and f f M (9, t) that of the amb ient ISM seen by 
stars of age t at galaxy inclin ation 9. Following IWild et"al] (2007, 
also lda Cunha et alj|2008h . we take the attenuation curve in the 
birth clouds to be 



U/0.55 nm)" 



(15) 



where f^ c is the angle- averaged V-band attenuation optical depth. 
We describe f' SM (#, i) using the multi-component model of 



T04. By analogy with Section [2731 we relate stars in different age 
ranges in a model galaxy to different geometric components of 
the T04 model. This is supported by observational evidence that 
stars in quiescent (i.e. non-starburst) star-forming galaxies tend to 
form along spiral arms in a thin disc close to the equatorial plane. 
The scale height of these stars can lat er increase as a consequence 
of dynamical heat ing and merging (jYoachim & Dalcantori 120061 ; 
Roska r et al]|2012h . as suggested by the distinct stellar populations 
of t he thin- and thick-disc components observed in ext ernal galax- 
ies dYoachim & Dalcan ton 2008; Yoachim et alj|2012n . Finally, in 
normal star-forming galaxies, the oldest s tellar populations are gen - 
erally found in a bulge component (e.g. IKauffmann et al. 2003b). 
We therefore consider the same three age parameters t,hi n , tthkk an d 



tbuige as in Section |2T3] such that t lhin ^ t^k ^ t bulge , and associate 
stars younger than t t hm with the thin-disc component, those in the 
age range tthi n < t ^ t t hi c k with the thick-disc component and those 
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D„4000 




D„4000 

Figure 16. (a) Balmer stellar absorption-line index HSa versus 4000-A dis- 
continuity D„4000 for the low-//, galaxies in the Wl 1 sample described in 
Section [3-1-11 (solid contours, encompassing areas of 1, 10 and 20 times 
the mean density) and the library of 5000 model spectral energy distribu- 
tions generated in Section 13.2.11 to analyse this sample (dots and dashed 
contours), (b) Same as (a), but for the high-//, galaxies. 



in the age range t t hk-k < t ^ tb u i g e with the bulge component. We 
compute the attenuation optical depth in the ambient ISM seen by 
stars of age / at galaxy inclination 9 in equation ( 114b as 



-ISM 

'a 



(,e,t) 



* bulge 



o 



for 
for 
for 
for 



t ^ tfhin , 
tthin < ^ ^ tthick , 
tthick < t ^ tbuige , 
t > tbuige > 



(16) 



where f* ln 



W), 



^thick^ 



(6) and f\ ulgl: (#) are the attenuation optical 



depths of the thin-disc, thick-disc and bulge components in the T04 
model. All three quantities are parametrized in terms of the central 
face-on B-band optical depth of the combined thin and thick dust 
discs, t b± (Section r2.1.2b . We note that the age parameters t t hm> ttMck 
and tbuige allow one to explore any combination of the T04 geomet- 
ric components: a model can be described by a single component 
(i.e. a pure thin disc, if tthin = tthick = tbuige; a pure thick disc, if there 
are no stars younger than t th i„ and if t thick = t b ui B e; and a pure bulge, 
if there are no stars younger than t th i„ = t thick ), two components (e.g. 
a thin and a thick stellar discs, if t th i n t t t hi C k and t, h ick = t bu i g e) or 
three components (if t thin + t thick # tbuige)- The precise choice of 
tthin! tthick and tb u ige, combined with the star formation history, de- 
termines the contributions by the different components to the total 
intensity produced by stars in the galaxy (weights ^' A ; equation 0. 
In summary, we model the attenuation of starlight by dust using 
five adjustable parameters, f^ c , t b± , t tMn , t thick and t bu i a e- 

In Table [2] we list the ranges in these parameters adopted in 
the next section to interpret observations of nearby star-forming 



galaxies. The largest observed Ha/Hfl ratio in the sample of Sec- 
tion [3J] corresponds to a total V-band attenuation optical depth 
of roughly 1.3. This includes both attenuation in the birth clouds 
of ionising stars and through the intervening ambient ISM. We 
therefore sample f^ c uniformly in the interval from to 1.5. Fur- 
thermore, we sample t b± uniformly across the full range consid- 
ered by T04, i.e. from to 8. Finally, to satisfy the condition 
tthin < tthick < tbuige, we take t thin (in Gyr) to be distributed log- 
arithmically between -3 and 1.13, t tn i c k logarithmically between 
max[log(t t hin/Gyr), —1.0] and 1.13, and tbuige logarithmically be- 
tween max[log(t t hick/Gyr), 0.0] and 1.13. 



3.3 A Bayesian approach: Markov Chain Monte Carlo with 
the Metropolis-Hastings algorithm 

We now use the library of model spectral energy distributions pre- 
sented in Section [3~2l to interpret the observed dependence of the 
Hcy/H/? ratio and ugrizYJH attenuation curve on geometry and in- 
clination in the sample of nearby star-forming galaxies described in 
Section uTTI To search for the combination of parameters fy C , t b± , 
tthin, tthick and tb u i 6e that best reproduces the observed trends, we 
adopt a Bayesian approach. Technically, this amounts to comput- 
ing posterior probability distributions of the adjustable parameters 
from the thorough exploration of the parameter space, given initial 
prior distributions and the data available to constraint the model. 
The Markov Chain Monte Carlo (MCMC) approach allows effi- 
cient exploration of complex, multi-dimensional parameter spaces. 
It differs from basic Monte Carlo methods through the inclusion 
of correlations between the random steps composing a chain, each 
step depending on the previous one (Markovian property). Such 
correlations allow an efficient exploration of the parameter space, 
making the algorithm spend most of the computation time in re- 
gions of highest probability (see Appendix [C]fbr more detail). 

To sample in a Markovian way the posterior probability dis- 
tributions of the five adjustable parameters of the dust model, we 
adopt the widely used random-walk Metr opolis-Hastings algorithm 
dMetropolis et al J 1 19531 : lHa stings 197 ^), as implemente d in the 
publicly available code CosmoMC jLewis & Bridle 2002). At each 
draw (step j) of a new set of parameters from a 'proposal' distribu- 
tion (see below), we compute the posterior distribution of this set 
by comparing the predicted and observed mean spectral properties 
of galaxies as a function of axis ratio. In practice, the adoption of 
uniform priors makes the posterior distribution proportional to the 
likelihood distribution, which we compute assuming independent 
Gaussian errors on the observations, as 



In this expression, the index k runs across the 10 bins of b/a and 
the index i across all spectral properties, from Ha/H/J ratio to dif- 
ferential ugrizYJH magnitudes. The quantity O*(0j) is the value 
of the i th observable obtained when using the set of parameters 
&'j in the dust model of Section [3,2,21 while O* is the observed 
(bias-corrected) value of that observable from Section 13.1.31 We 
estimate the observational error <r* from the variance of the im- 
portance weights used to correct the data for biases (Section [3. 1.31 
and Appendix|Bj we add an additional error of 0.02 mag in quadra- 
ture to account for the residual biases for the differential ugrizYJH 
magnitudes; such a correction is negligible for the Hfy/H/? ratio). 
The ratio of the posterior probabilities of two consecutive sets of 
parameters orients the move toward the high-probability region of 
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Parameter Definition Lower bound Upper bound 



'v 


V-band attenuation optical depth of dust in 
stellar birth clouds. 


0.0 


1.5 


TB± 


Central face-on B-band attenuation optical 
depth of dust in the diffuse ISM. 





8 


log(t thin /Gyr) 


Stars with ages t ^ t tn jn are associated with 
the T04 thin-disc component. 


-3.0 


1.13 


log(t t hick/Gyr) 


Stars with ages t t hin < f ^ 'thick are associ- 
ated with the T04 thick-disc component. 


max[log(t thin /Gyr), - 1 .0] 


1.13 


l°g(tbui ge /Gyr) 


Stars with ages tjhici < t ^ tb u i„ e are associ- 
ated with the T04 bulge component. 


maxIlogCtthiek/GyrJ.O.O] 


1.13 



Table 2. Prior distributions of the five main adjustable parameters of the dust model described in Section l3.2.2l 



the parameter space. The result of the algorithm is a sequence, or 
chain, of sets of parameters 0^, which build up step by step the 
posterior probability density functions of f^ c , Tb±, tthin. tthick an d 

tbulge- 

A more detailed description of this procedure is presented in 
Fig. [T7] The chain starts at step j = (not shown) with a random 
drawing of parameters within the ranges given in Table [2] At each 
new step j, a candidate set of parameters 0^ is drawn from the 
proposal distribution, which in CosmoMC is a multi-variate nor- 
mal distribution, noted AA(0 7 -i,C), centred on the previous set 
7 _i and with covariance matrix crj Then, we use the candidate 
set of parameters 0^ to attenuate the spectral energy distributions 
of the 5000 stellar populations in the library assembled in Sec- 
tion [X2J] as described in Section [3.2.21 This allows us to com- 
pute the likelihood -£(0^) of the candidate set by comparing the 
predicted and observed mean spectral properties of galaxies as a 
function of b/a (we assume Gaussian observational errors). The 
candidate set is accepted, i.e. we adopt ®j = 0' ; , if the likelihood 
ratio X(0j)/X(0j-t) is larger than a critical value r drawn uni- 
formly at random from between and 1. Otherwise, we draw a 
new candidate set of parameters 0^. 

We perform these steps along six chains run in parallel, un- 
til stationary posterior probability density functions are reached 
(this is usually achieved for j around a few thousand). We de- 
termine convergence through a comparison of the intra- and 
inter-chain variance for each parameter, setting the criterion of 
iGelman & Rubinl Jl992h to 1 - R = 0.01. As is standard with 
CosmoMC, we compute final posterior probability distributions by 
combining the samples produced in the last half of all chains. In 
Appendix IC2I we provide more detail on this procedure, as well 
as a visual representation of the build-up and adequate sampling of 
the posterior probability distribution of each parameter. 



Pick new set of parameters 

0) «- N(&j-i,C) 



Proposal 
covariance 

c 



1 



1 



Diffuse ISM 



Birth clouds 



Total attenuation 
.to 



tf(0,t) 



Attenuation of library of 
galaxies spectra 



Compute likelihood 



, if vfte'j) 



Chariot & Brazual 
templates 



Star formation and 
chemical enrichment 
histories 



Draw random number 
r e [0, 1] 

£(0'j) 
if > r then 



} = 0) 

else j - 0j_\ 



Figure 17 . Flow chart summarising the random-w alk Metropolis-Hastings 
algorithm ^Metropolis et al. 1953; Hastinas 1970) adopted in Section [3~4l to 
sample the posterior probability distributions of the adjustable parameters 
^y C > r s±, tthin, tthick and tb u ige, assuming flat, truncated priors in the ranges 
given in Tablef2] 



Theoretically, the optimum choice for C is some scaled ver- 
sion of the covariance ma trix of the posterior probability distributions 
(Roberts & Rosenthal 200l]). In practice, this quantity being unknown in 
our case, the possibility with CosmoMC to build up C th rough step- 
by-ste p adjustments is extremely time consuming. Following iTrotta et alj 
201 lh. we therefore appeal to the fre ely available code MultiNest 
Feroz & Hobson 2008; Feroz et alj|2009l) to obtain a fast estim ate of C 
based on the nested sampling algorithm of Skilling (2004, 2006, see Ap- 
pendix|D]for detail). 



3.4 Constraints on the content and spatial distribution of 
dust in nearby star-forming galaxies 

Fig. ED shows the posterior probability distributions of t^ c , t b ±, 
tthin, tthick and tbnigc obtained from the MCMC analysis of the Hcf/H/J 
ratios and ugrizYJH attenuation curves of the sample of 1ow-/j,, 
nearby star-forming galaxies of Section|3T| Table|3]lists the corre- 
sponding median and 2.5th-97.5th percentile range of each param- 
eter. To complement this information, we compare in Fig. |19|the 
mean Hcr/H/J ratio and differential ugrizYJH magnitudes predicted 
using the best-estimate (i.e. median) parameters from Fig.[T8] with 
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0.0 0.5 
log(t,hin/Gyr) 



0.4 0.! 
log(t th ick/Gyr) 



0.0 0.5 1.0 
log(t bu ige/Gyr) 



Figure 18. Results of the Bayesian MCMC analysis of the dependence of the Hct/H/? ratio and ugrizYJH attenuation curve on galaxy inclination in the Wl 1 
sample of nearby star-forming galaxies described in Section 13.1.11 using the library of model spectral energy distributions assembled in Section 13.21 The 
diagonal panels show the marginal posterior probability distributions for f^ c , rg ± . t t hm> t t hick and tbuige, for the low-/i, (blue lines) and high-/j, (red lines) 
galaxies. The off-diagonal panels show the joint probability distributions of the same parameters, for the low-/;, (blue contours) and high-/i* (red contours) 
galaxies. In each of these panels, the contours enclose 68 and 95 percent of the posterior probability. 



the (bias-corrected) observations of Fig. QT] as function of galaxy 
inclination. The agreement between model and data is remarkable 
at all wavelengths in this figure. 

We now describe these results for the 1ow-/j„ galaxies. The 
blue lines on the diagonal panels of Fig. Q~8] show the marginal 
posterior probability distributions of the five dust parameters. For 
f^ c , t b± , t thin and t thick , the distributions are nearly Gaussian with 
well-defined peaks and tails. For t bu i gc , the distribution peaks at 
log(t bll i„ c /Gyr) = 1.1 and is truncated at the maximum allowed 
age of log(t bu i gc /Gyr) = 1.13 (t bulge = 13.5 Gyr). To interpret the 
marginal distributions of Fig. II 81 (Table l3t. it is helpful to return 
to the definition of the age parameters in Table [2] In this context, 
we infer from the median values in Table[3]that stars younger than 
t as 2.6 Gyr [log(//Gyr) as 0.42] in nearby low-/(» galaxies are atten- 
uated in the same way as T04 thin-disc stars, those in the age range 
2.6 < t < 3.9 Gyr [0.42 < logO/Gyr) < 0.59] in the same way as 
T04 thick-disc stars and those in the age range 3.9 ^ t < 12.6 Gyr 
[0.59 $J log(//Gyr) < 1.1] in the same was as T04 bulge stars. It 
is also useful to remember that tthin, tthick and tbuige* when combined 
with the star formation histories of model galaxies (Section |3. 2b . 



determine the intensity weights f ', of the thin-disc, thick-disc and 
bulge components in the attenuation arising form the diffuse ISM 
(see equation[3}- Hence, we can relate t t hi n , tthick and tbuige to V-band 
intensity weights £* ick and ^ ulge , which we report in Table [3] 
The results suggest that 80 percent of the V-band attenuation in 
the diffuse ISM of low-//, galaxies in our sample is characteris- 
tic of that affecting thin-disc stars in the T04 model, 9 percent to 
that affecting thick-disc stars, and the remaining 10 percent to that 
affecting bulge stars. We note that, according to the results of Sec- 
tion |3.6| below, the finding of a non-zero contribution by bulge stars 
to the integrated light of low-//, galaxies, which are expected to 
be bulge-less (Section [3.1. It , could arise from an imperfect cor- 
rection of systematic biases in the dependence of dust attenuation 
signatures on inclination in the Wll sample (Section |3. 1 .3i . The 
median t b± » 1.9 favoured by our analysis of low-//, galaxies 
corresponds to an angle-averaged V-band attenuation optical depth 
("fy }e ~ 0.26 in the ambient ISM, and the median attenuation op- 
tical depth of stellar birth clouds is f^ c as 0.48. 

The bottom-left off-diagonal panels of Fig. [T8l show the joint 
posterior probability distributions of all pair of parameters. The 
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joint probability distribution of f ^ c and t b± does not show any 
degeneracy, indicating that the combined analysis of the Hcr/H/? 
ratio and ugrizYJH attenuation curve allows us to efficiently sep- 
arate the attenuation arising from the diffuse ISM from that per- 
taining to stellar birth clouds. The slight apparent degeneracies 
affecting the age parameters t thick and t thin and the tails of t blllgc 
and t t hi c k, arise from the presence of likely models at the limit of 
the tthta < tthidj ^ tbuige criterion (i.e. around t thm = t thick and 
tthick = tbuig C ; see prior definitions in Table|2}. 

Figs[T8]and|20]show the analog of Figs[T8]and[T9]for the high- 
/j» galaxy sample. As before, we report in Table [3]fhe median and 
2.5th-97.5th percentile range of the marginal posterior probability 
distribution of each parameter. In Fig. [20] the agreement between 
predicted and observed mean spectral properties as a function of 
galaxy inclination is again remarkable. 

As in the case of low-/i, galaxies, the marginal distributions 
of f^ c , t b ±, tthin, tthick and tb u ig C shown as red lines in Fig. [18] are 
nearly Gaussian with well-defined peaks and tails and with a cut- 
off for thuige at the maximum age log(t bu ig C /Gyr) = 1.13. The me- 
dian parameter values favoured by observations of high-yu» galax- 
ies in Table [3] differ slightly from those found above for low-/i, 
galaxies. The slightly larger t, h i n (4.7 versus 2.6 Gyr), t thick (8.2 ver- 
sus 3.9 Gyr) and smaller tb u i ge (8.7 versus 12.6 Gyr) imply small 
changes in the K-band intensity weights £*"*and ^ ulge re- 
ported in Table[3] The results are similar to those for low-yu, galax- 
ies in that the V-band attenuation of the diffuse ISM of high-yu„ 
galaxies is primarily characteristic of that affecting thin-disc stars 
in the T04 model, while 16 percent is more similar to that affecting 
thick-disc stars, and 3 percent bulge stars. Also, we find a median 
t b ± ~ 1.7, corresponding to an angle-averaged V-band attenuation 
optical depth (t£ m }# ~ 0.31 in the ambient ISM, very similar to the 
value found for low-yU» galaxies. In contrast, the favoured median 
attenuation optical depth of stellar birth clouds is roughly twice 
larger for high- than for low-yU» galaxies (f^ c « 0.85 versus 0.48). 

The top-right off-diagonal panels of Fig. [TS] show the joint 
posterior probability distributions of all parameters, revealing a 
weak degeneracy between f^ c and t b± . As for the low-p, galaxies, 
the joint posterior probability distributions of the age parameters 
tthin, tthk k and tbuigc appear degenerate because of the presence of 
likely models close to the t t hhi ^ Wk ^ tbuige limit. 

We have now established that the model presented in Section 
13.21 can account remarkably well for the observed dependence of 
the Hff/H/3 ratio and ugrizYJH attenuation curve of both low- and 
high-/j, galaxies on inclination (Figs [T9land|20l. It is worth recall- 
ing that the predictions of the T04 model adopted in our MCMC 
analysis are similar to those of other types of state-of-the-art radia- 
tive transfer models (Section[272]and[23). Therefore, we do not ex- 
pect the results presented in this section to depend critically on the 
details of the radiative transfer model. These results have several 
important implications for our ability to constrain the content and 
spatial distribution of dust in nearby star-forming galaxies. Firstly, 
the ability with a single best estimate of the central face-on fi-band 
optical depth t b± to reproduce the dependence of the Ha/Hfl ratio 
and ugrizYJH attenuation curve on inclination across the full ob- 
served range in b/a implies that geometric effects alone can account 
for these observational trends. Also, the similarity of the best esti- 
mates of r B± for low- and high-p, galaxies (1.9 versus 1.7) implies 
that the surface mass density (and hence the mass, since both sam- 
ples have similar size distributions; see Fig. 1100 of diffuse-ISM 
dust should be roughly similar in both samples. 

Another important implication of the agreement between 
model and data in Figs [19] and [20] is that it validates, a posteri- 
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Figure 19. Dependence of integrated spectral properties on galaxy inclina- 
tion, as predicted by models with f^ c , t b± , t t hin, tthick and t bu i 6e parameters 
corresponding to the median values reported in Tablef3]for low-/j, galaxies 
(black lines), compared with the observed (bias-corrected) dependence in 
the Wll sample from Fig. II II (open circles). The grey shading locates all 
models with parameters spanning the 2.5th-97.5th percentile ranges of the 
posterior probability distributions of TablefJ|for low-//* galaxies. 



ori, our combination of the T04 model with a galaxy spectral evo- 
lution code by associating stars in different age ranges with the 
thin-disc, thick-disc and bulge components of T04. This demon- 
strates the possibility to relate integrated spectral properties of 
structurally unresolved galaxies, such as the Ha/H/3 ratio and dif- 
ferential ugrizYJH magnitudes, to different geometric components 
and the spatial distribution of dust in these galaxies. For example, 
we find here that over 80 percent of the K-band attenuation in the 
diffuse ISM of both low- and high-yu„ star-forming galaxies is char- 
acteristic of that affecting thin-disc stars in the T04 model. This by 
itself means that the luminosity of these galaxies is dominated by a 
component in which the spatial distribution of stars and dust have 
the same scale heights, as in the thin-disc model of T04. Finally, the 
fact that the best estimate of f ^ c is about twice as large for high- 
than low-fi, galaxies, while that of t b± is similar, sets interesting 
constraints on the distribution of dust in the ISM. This could sug- 
gest that, for example, low-yu, galaxies, which are typically more 
metal-poor and hence have presumably lower dust-to-gas ratio than 
their high-/j, counterparts (Fig. UOt , contain proportionally more 
gas than high-yu, galaxies, while the masses of individual stellar 
birth clouds are more similar. Since both types of galaxies have 
similar specific star formation rates (Fig. 1101 . this would also re- 
quire that low-yu, galaxies form stars less efficiently than high-yu, 
galaxies out of the available gas. 



3.5 Integrated versus fibre properties 

The constraints derived above on the content and spatial distribu- 
tion of dust in nearby star-forming galaxies pertain to an inner re- 
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Original model parameter Associated physical parameter Definition 

-BC low-//, 0.48+q J52 , BC 1ow-/j, °- 48 !ao2 V-band attenuation optical depth of dust in stellar birth 

Ty high-//, 0.85+002 T v high-//, 0.85+jJ-^ c i ou d s . 

low-//, l-9°-oju r isM\ low-//, 0.26+0^ Angle-averaged V-band attenuation optical depth of 

TS± high-//, 1.68+0 |» V high-/,, 0.31!0 04 dust in the diffuse ISM. 

, low-//, 0-42+0;|] low-//, 0.80+0;°^ y_b an d intensity weight of thin-disc stars in the T04 

logcwtxyrj h . gh _^ t 06? ^ ? v 0.80+°,^ model. 

. . low-//, 0.59+0 .18 ^ low-//, 0.09^0;^ y.band intensity weight of thick-disc stars in the T04 

log(t thlck /Uyr) h . gh _^ o9i+Q ^ f y oi6 ^^ modd 

, low-//, 1-08+0;^ ^ low-//, O.IO+O;] 2 , y_b an d intensity weight of bulge stars in the T04 

10g(t bulge /uyr, h . gh _^ Q94 ^^ ^ h . gh _^ oQ3 ^ modd 



Table 3. Median and 2.5th-97.5th percentile range of the posterior probability distributions of the five adjustable model parameters of Table [2] derived from 
the MCMC analysis of the low- and high-/;, galaxies in the Wll sample described in Section [TTTT] using the library of model spectral energy distributions 
assembled in Section [T2l 
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Figure 20. Same as Fig. 1191 but for the high-//, galaxies. The (bias- 
corrected) observations are from Fig. 1121 



gion of 3-arcsec diameter probed by the SDSS fibre. To check the 
applicability of these results to entire galaxies, we first build inte- 
grated ugrizYJH spectral energy energy distributions of the galax- 
ies in the Wll sample by computing, for each waveband, the flux 
comprised in a circular aperture of size the r-band Petrosian ra- 
dius (footnote |6). In practice, we interpolate the flux within this 
radius from observations in different apertures (available at 8 radii 
between 0.23 and 1 1.42 arcsec for SDSS data, and 12 radii between 
0.5 and 10.0 arcsec for UKIDSS data; we exclude ~ 5 percent of 
galaxies with r-band Petrosian radius larger than 10 arcsec in the 
original sample; see Fig.llOU). Then, we follow the procedure out- 
lined in Section [3.1.1l and compute the dependence of the attenua- 
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Figure 21. Difference between the mean magnitude of galaxies in bins of 
different inclination, from face-on (l-b/a » 0) to edge-on (l-b/a » 1), and 
that of face-on galaxies, at fixed K band flux (i.e. fixed stellar mass), for the 
low-//, (top) and high-//, (bottom) galaxies of the Wll sample described in 
Section l3.1.1l in the ugrizYJH bandpasses, as indicated. Different symbols 
refer to different aperture sizes. Open diamonds refer to differential mag- 
nitudes measured in a large circular aperture of size the r-band Petrosian 
radius. Filled circles with error bars refer to differential magnitudes mea- 
sured in a 3-arcsec-diameter aperture and corrected for the bias in galaxy 
size Rgo using the procedure described in Section |3.1.3| and Appendix |T?1 



tion curve on galaxy orientation, by comparing the mean total mag- 
nitudes of galaxies in different bins of b/a to those in the face-one 
bin. 

We do not have access to spectroscopic measurements out- 
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side the SDSS-fibre aperture. Thus, unlike in Section l3.1.1l we can 
neither measure the dependence of the integrated Ha/H/j ratio on 
inclination, nor correct the dependence on inclination of the total 
ugrizYJH attenuation curve for systematic biases in gas-phase oxy- 
gen abundance and specific star formation rate. Yet, we can com- 
pare in a meaningful way fibre versus total dust properties by cor- 
recting the dependence on inclination of the fibre-derived attenua- 
tion curves purely for the bias in galaxy size Rgo, and not applying 
any correction to the total attenuation curves. This is because the 
fraction of total light sampled by the SDSS fibre depends on galaxy 
size, which in the Wll sample correlates with axis ratio (Figs 1 131 
and 1 14b. while magnitudes measured within the r-band Petrosian 
radius should include the majority of the light of a galaxy, at any 
inclination. 

Fig. [21] shows the differential ugrizYJH magnitudes, Am = 
m(b/a) - m(face-on), obtained in this way as a function of b/a, for 
the low- and high-p, samples, respectively. Open diamonds refer 
to the total dust attenuation curves, computed directly from obser- 
vations. Filled circles refer to the fibre-derived attenuation curves, 
corrected for the bias in galaxy size Rgo using the procedure de- 
scribed in Section [3. 1.3| and Appendix |B1 The remarkable agree- 
ment, within the errors, between the total and fibre-derived attenua- 
tion curves in Fig.|2T|suggests that the results derived in Section [3~4l 
(Table [3} from fibre observations are applicable to entire galaxies. 
We note that this agreement in differential ugrizYJH magnitudes 
cannot reveal potential offsets in the absolute normalisations of the 
attenuation curves. In fact, Wll find that the V-band attenuation 
optical depth of galaxies in their sample is typically ~ 15 percent 
larger in the central regions than in an aperture of size R 90 , indepen- 
dent of axis ratio and stellar-mass density (see their fig. 6). Hence, 
the constraints on (f^ M )s derived from our analysis of fibre obser- 
vations in Section [3~l4l may slightly overestimate the galaxy-wide 
attenuation of Wll galaxies. 

3.6 Impact of the correction of systematic biases in the 

dependence of dust attenuation signatures on inclination 

The results reported in Figs ll84j2Ql were obtained after correcting 
the observed dependence of the Ha/H/? ratio and ugrizYJH atten- 
uation curve of nearby star-forming galaxies on inclination for var- 
ious systematic biases, as described in Section [3. 1.3 1 To quantify 
the impact of these corrections on our results, we have repeated the 
same analysis as in Section [3~l4l but fitting this time the uncorrected 
data (open circles in Figsll llandll2t in place of the corrected ones 
(filled circles). 

We find that, in this case, the model presented in Section [3T21 
can still account reasonably well the observations of 1ow-/j, galax- 
ies in Fig. QT| although the best-estimate (i.e. median) parame- 
ters differ from those derived in Section |3~4] The favoured optical 
depths are r s± as 0.41 and ty C ss 0.59, and the age parameters 
tthin ~ 0.02 Gyr, t thick «0.32 Gyr and t bu i gc « 8.0 Gyr indicate that 
the T04 bulge accounts for 60 percent of the V-band attenuation in 
the diffuse ISM. The reason for this is that the ugrizYJH attenua- 
tion curve of uncorrected data depends only weakly on galaxy ori- 
entation in Fig.Qj] which is characteristic of the T04 bulge model. 

Similarly, the model can still account well for the observations 
of high-yu„ galaxies if no bias correction is applied to the data. The 
best-estimate optical depths, t b ± ~ 2.3 and f^ c ss 0.78, are close 
to the values found above when accounting for bias corrections. 
The constraints on the age parameters are slightly different, corre- 
sponding to t t hi n ~ 3.0 Gyr, t th kk ~ 4.4 Gyr and t bl ,i gc a; 4.8 Gyr. The 
fact that all three ages are close together indicates that the V-band 



attenuation in the diffuse ISM is almost exclusively characteristic 
of that affecting thin-disc stars in the T04 model. This makes the 
fraction of V-band luminosity produced by unattenuated stars, with 
ages greater than tbuige. increase from roughly to 20 percent. 

Hence, as expected from the stronger bias corrections derived 
in Section [3.1.3l for low- than high-yu, galaxies, the changes in pa- 
rameter estimates induced by these corrections are significantly 
larger for low- than for high-yu, galaxies. 



4 DISCUSSION AND CONCLUSION 

4.1 General results on dust attenuation curves from models 
of radiative transfer in dusty media 

We have presented a new approach to investigate the content and 
spatial distribution of dust in structurally unresolved star-forming 
galaxies from the dependence of integrated spectral properties on 
inclination. To achieve this, we first established that four differ- 
ent types of existing models of radiative transfer in disc galax- 
ies all predict a similar dependence of the dust attenuation opti- 
cal depth in the ambient ISM, f' SM (#), on viewing angle 9 in the 
ugrizYJH bandpasses. We find that, in fact, these models predict a 
quasi-universal relation between slope of the attenuation curve at 
any wavelength from the ultraviolet to the near infrared, nf M (G), 
and V-band attenuation optical depth of the dust, f¥? M (0), which we 
parametrize by means of a simple analytic expression (equation [9] 
of Section [2j2j. As shown by Figs|4^ and[4ji, the range in n' SM (#) 
spanned by this relation at A = 0.55 and 1 .6 jim acro ss the ful l range 
in f'y M (0) is much wider than that sampled by the lCalzettil d200ll) 
attenuation curve and the Milky-Way, LMC and SMC extinction 
curves, even though all the models considered in the present study 
rely on optical properties of Milky Way-type dust. This confirms 
that, in normal star-forming (disc) galaxies, changes in the shape 
of the optical attenuation curve induced by geometry and orienta- 
tion effects are likely to dominate over those in duced by differences 
in the optical prope rties of dust grains (e.g., iGranato et alUkOOCt 
iFontanot et al]|2009h . 

The quasi-universal relation between n' SM (#) and t l y M (6) can 
be easily implemented in any spectral analysis of structurally un- 
resolved disc galaxies to account in a simple yet physically con- 
sistent w ay for the effect of galaxy orientation on dust attenuation 
(see e.g. |Pacificietai]|2012h . This relation may also be relevant 
to the interpretation of recent observation al results by Wuvts et ail 
(2011) and Arnouts et. al (in preparation). IWuvts et al.l d201 lb find 
that, at redshifts z < 3, the rest-frame ultraviolet, optical and 
near-infrared col ours of actively star-forming galaxies can be best 
modelled using a lCalzettil J200lh attenuation curve, while those of 
more quiescent star-formin g galaxies are bet ter reproduced with a 
steeper, SMC-like curve. IWuvts et alj d201ll) suggest that system- 
atic changes in the grain size distribution could account for this 
effect. A similar correlation between shape of the optical attenua- 
tion curve and star formation activity has been found by Arnouts et. 
al (in preparation) in a sample of z < 1.5 star-forming galaxies. We 
argue here that behaviours of this type could result naturally from 
the relation between nf M (9) and Pyiff), if star formation activity 
correlates with dust attenuation optical depth in the ambient ISM 
of galaxies. Such a correlation is expected from the combination 
of the Schmidt-Kennicutt relation jSchmidt|[T95^ ; lKennicutt1l 19981) 
with the relat ion between metallicit y and star formation rate at fixed 
stellar mass (Ma nnucci et alj|2010l) . According to these trends, ac- 
tively star-forming galaxies should contain more gas and metals, 
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and hence presum ably more dust tha n more quiescent star-forming 
galaxies (see also lZahid et alj|2013h . It is possible that, therefore, 
the observed systematic changes of galaxy attenuation curves with 
star formation activity might arise purely from geometric effect, 
rather than from changes in the physical properties of dust grains. 

4.2 Beyond angle-averaged dust prescriptions 

These findings demonstrate the importance of accounting for the 
dependence of integrated spectral properties on inclination when 
constraining galaxy physical parameters. In this paper, we have 
proposed an original, computationally inexpensive approach to ac- 
complish this goal by extending simple angle-averaged models of 
dust attenuation using the generic output of sophisticated radia- 
tive transfer calculations . Specifically, we have extended the simple 
angle-averaged model o flCharlot & Faiil feoOO), which accounts for 
the different attenuation affecting young stars in their birth clouds 
and older stars in the ambient ISM, by associating stars in different 
age ranges in the ambient ISM with different geometric compo- 
nents of galaxies (Section [3.2.21 . This is motivated by the obser- 
vation that different stellar populations reside in different spatial 
components of nearby star-forming galaxies. Stars in these galax- 
ies tend to form along spiral arms in a thin disc close to the equa- 
torial plane (as a consequence of gas dissipation). The scale height 
of these stars can later increase as a result of dynamica l heating 
and merging ( Yoachim & Dalcanton 2006; Rosk ar et al.ll2012h . as 
suggested by the presence of different ste llar populations in the 
thin and thick discs of ex ternal galaxies dYoachim & Dalcanton! 
120081 : lYoachim et alj|2012|). while th e oldest stars reside typically 
in a bulge iKauffmann et alj|2003bl) . The association of different 
stellar-age ranges with different geometric components of galaxies 
presents several notable advantages with respect to traditional in- 
vestigations of spatially resolv ed nearby galaxies using detailed ra- 
diative transfer modelling (e . gjBaes et al.l20 10; MacLachla n et alj 
l201ll ; ldeLoozeet"afll2012h |Pe Gevter et alj|2013l) . Firstly, it can 
be easily combined with sophisticated models of spectral evolu- 
tion. Also, it is applicable to the study of structurally unresolved 
galaxies. And finally, it can be efficiently applied to the statistical 
analysis of large samples of galaxies. 

We have implemented this approac h by combi ning the multi- 
component radiative transfer model of iTuffs et alj ^20041) with a 
comprehensiv e library of star form ation and chemical enrichment 
histories fromlPacifici et al.|( 2012 ) and a stellar population synthe- 
sis code ( iBruzual & Ch ariot 2003; Chariot & Bruzual, in prepara- 
tion). The resulting library of galaxy spectral energy distributions 
allowed us to investigate the dependence of dust attenuation on in- 
cli nation in the sampl e of nearby star-forming galaxies assembled 
bv lWild etaLU20T lb). To avoid unwanted biases in the observa- 
tional measurements, we also developed an original method, based 
on importance sampling, to measure the dependence of the Hcr/H/? 
ratio and ugrizYJH attenuation curve as a function of galaxy axis 
ratio in this sample (Section[33}- From the Bayesian MCMC anal- 
ysis of two subsamples of galaxies with low (bulge-less) and high 
(with a bulge) stellar-mass density ji,, we can draw valuable con- 
clusions on the content and spatial distribution of dust in nearby 
star-forming galaxies, which are encoded in the posterior probabil- 
ity distributions of the five main adjustable model parameters: the 
V-band attenuation optical depth of stellar birth clouds, f ^ c (as- 
sumed independent of angle); the central face-on B-band optical 
depth of dust in the ambient ISM, t b± , and the age parameters t thin , 
tihick an d tbuige separating stars in the thin-disc, thick-disc and bulge 
components (Tables [2] and [3). Our results suggest that fy C is sig- 



nificantly larger in high-yu, than in low-yU„ galaxies, while r B± is 
roughly similar in both cases. As mentioned in Section 13.41 this 
can arise if, for example, high-yu„ galaxies have higher star forma- 
tion efficiency than their low-yu, counterparts at fixed specific star 
formation rate, and dustier stellar birth clouds because of the higher 
metallicity (Fig.llOll. 



4.3 Content and spatial distribution of dust in the Wll 
sample of star forming galaxies 

The central face-on S-band optical depth of dust in the ambient 
ISM favoured by our analysis of nearby star-forming galaxies is 
T fl± « 1.8*2'? (Tabled- This is slightly larger than the range 0.5 < 
tb± ^$11 derived from detailed modelling of the BVRIK vertical 
surface brightness p rofiles of seven edge-on late-type galaxies by 
IXilouris et all fl999). At the same time, the value favoured by our 
analysis is smalle r than the constraint t b± ss 3.7 ± 0.7 inferred by 
iDriver et alj ( 120071) from an analysis of the local B-band luminosity 
function of galactic discs. These authors measure the dependence 
of the turnover of the galaxy luminosity function on inclination and 
interpret the fading of the turnover luminosity fr om low to high in- 
clination as a signature of attenuation by dust. IDriver et alj ( 2007) 
find that the observations can be well reproduced by a T04 thin- 
disc model with r B± =s 3.7 (we note that, according to their fig. 10, 
a value r B± =s 2 appears to also provide excellent agreement with 
their data). An advantage of the approach developed in the present 
paper over this study is that it makes use of the entire spectral in- 
formation in the ugrizYJH bandpasses. 

Our results on dust attenuation in nearby star-forming galax- 
ies are also consistent with the conclusions from recent studies 
of low-redshift galaxies observed at infrared wavelengths with the 
Herschel Space Observatory. In particular, based on the spectral 
analysis of ultraviolet, optical and infrared spectral energy distri- 
bution s of a sample of about 300 nearby galaxies, ICortese et alj 
j2012l) find a strong anti-correlation between dust-to-stellar mass 
ratio, Md llst /M,,, and stellar surface mass density, ji,. This relation 
supports our finding that low-yu, galaxies, which have lower stel- 
lar masses than high-/i» galaxies (Fig. [lOb). also have similar dust 
masses (as probed by the central face-on 6-band optical depth r B± ; 
Section [3~4l . We note that the similar distributions of specific star 
formation rates of low- and high-yu» galaxies in the Wll sample 
(Fig.llOb) prevent u s from testing the correl atio n between M d ust/M, 
and if/,; identified bv lda Cunha et alj 1 I2OIOI) and lSmith etal] j2012h 
in two independent samples of about 1500 nearby star-forming 
galaxies. Our results are consistent with the large scatter about this 
mean relation. 

We also note th at, during the revision of this manuscript, 
iGrootes et all d2013l) proposed a relation between t b± and 
stellar surface mass density ji,, of the form logr fl± = 
1.12(±0.11) log(/i,/M pc- 2 ) - 8.6(±0.8), based on the analysis 
of the ultraviolet, optical and infrared spectral energy distribu- 
tions of a small sample of 74 nearby disc galaxies. This implies 
r B± = 3-0*11 for low-yU„ galaxies in the Wll sample [median 



log(//,/MQpc~ 2 ) = 8.1] and t b± = 13.9^ 2 for high-^/, galaxies 
[median log^./MQpc -2 ) = 8.7]. The value for high-^/, galaxies is 
much larger than that inferred in our analysis (and in any of the pre- 
viously quoted studies of nearby disc galaxies). The origin of this 
discrepancy is unclear. It may arise, at least in part, f rom the differ- 
ent ran ge in i/r s spanned at fixed /i, by galaxies in the lGrootes et alj 
j2013l) sample compared to the Wll sample. We note that the dis- 
persion in r B± (t b in their notation) at fixed ji, exceeds an order 
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of magnitude in the small sample of nearby galaxies studied by 
iGrootes et all d201 3l) . 

In our approach, the age parameters t thin , t thick and t bll i ec , along 
with the adopted library of star formation histories, determine the 
relative weights of the different geometric components of galaxies 
in the overall dust atten uation curve. We find that, for both low- and 
high-yu, galaxies in the lWild et al] d2011bh sample, over 80 percent 
of the attenuation is characteristic of that affecting thin-disc stars 
in the T04 model, while the remaining 20 percent is more similar 
to that affecting t hick-disc and bulge stars (Table |3j we note that 
IWild et al]|2011bl also suggest that their observations might be ex- 
plained by a scale height for the dust that is smaller than that of 
old stars). It is interesting to examine the ratio of thick-disc to thin- 
disc stellar luminosity predicted by our analysis, since this quantity 
has been estimated observationally. Our results suggest that the me- 
dian unattenuated V-band luminosity ratio of thick-disc to thin-disc 
stars is fS" /£y ~ 0.1 an d 0.2 for low- and high-^u, ga laxies, re- 
spectively. Observationally. lYoachim & Dalcantor] ( 120061) find that 
the i?-band luminosity ratio of thick-disc to thin-disc stars (cor- 
rected for internal extinction) in a sample of 34 nearby edge-on disc 
galaxies decreases as a function of circular (i.e. disc rotation) ve- 



locity roughly as tf^/i 


-* in = 0.25(K/100kms-')- 21 (theirfig. 18; 


see alsolComeron et al. 


2012b. We can estimate the circular veloc- 


ities of galaxies in the 


Wild et alj J2011bh sample from the avail- 



J201 ]]). We find median circular velocities V c « 150 km s 1 and 
190 km s -1 for low- and high-p, galaxies, respectively (correspond- 
ing to median masses 6.3 x 10 9 Mq and 3.2 x 10 10 Mq), which im- 
plies ^j 1 "*/^*" 1 ~ 0.1 for both types of galaxies. This is in remark- 
able agreement with the results of our analysis, given the negligible 
differences expected between V- and W-band luminosity ratios. 



4.4 Effect of galaxy inclination and spatial distribution of 
dust on broadband colours 

The much wider range of slopes of attenuation curves probed by the 
models in Fig. [4] relative to commonly adopt ed dust prescr iptions 
(e.g. Milky- Way and SMC extinc tion curves; ICalzettill200ll atten- 
uation curve; ICharlot & Falll2000l two-component dust model) has 
important implications for galaxy spectral analyses. We can illus- 
trate this by comparing the spectral energy distributions of galax- 
ies computed using standard dust prescriptions with those derived 
using equations f7t-dl It above, which account for the effects of 
galaxy inclination on the dust attenuation curve. For simplicity, we 
appeal to the same library of 5000 star formation and chemical en- 
richment histories as used in Section [3j2JJ ab ove to describe the 
properties of 1ow-/j, galaxies in the lWild et al.1 ( feOllbl) sample. We 
compute attenuation by dust in these galaxies using different pre- 
scriptions. For each galaxy in the library, we first draw the V-band 
attenuation optical depth in the ambient ISM from a Gaussian dis- 
tribution centred on f'^ M = 0.65 with a dispersion of 0.40. Then, 
we draw the parameters »y M and b defining the attenuation curve 
(see equation [9]l from Gaussian distributions centred on the val- 
ues predicted by equations d 1 0| > and ( lilt , with dispersions of 25 
and 10 percent, respectively. Fig. |22]shows the difference between 
the ugrizYJHK colours of galax ies computed u sing this prescrip- 
tion and those obtained using the lCalzettil d200lh attenuation curve 
corresponding to the same V-band attenuation optical depth, for 
each galaxy (black points). The differences are substantial, reaching 
0.3-0.4 mag in optical colours and 0.1 mag in near-infrared ones. 
To compare our prescr iption with the two-component dust model 
of lCharlot&Falll ( l2000l) , we also include attenuation by dust in 
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Figure 22. Difference in the broadband colours of 5000 model galaxies 
computed using the same star formation and chemical enrichment histo- 
ries as in Fig. l 16t . but different prescriptions for attenuation by dust. Black 
points: colours computed using the approach presented in Section 13.2.21 
above, which account for the effects of galaxy inclination on the dust at- 
tenuation curve, minus colours computed with the prescriptions of lCalzettl 
(2001). Grey points: colours computed using the approach presented in Sec- 
tion sec:d ust above, minus colou rs computed using the two-component dust 
model of Chariot & Fall 12000). In both cases, the substantial differences 
illustrate the importance of accounting for the influence of geometry and 
orientation effects on galaxy spectral energy distributions. 



stellar birth clouds by drawing the ratio f!f M /(f!f M + f^ c ) from a 
Gaussian distribution centred on 0.4 with a dispersion of 0.2. The 
colour differences in Fig.|22]in this case are again substantial (grey 
points). The results of these tests indicate that neglecting the ef- 
fect of geometry and orientation on attenuation can severely bias 
the interpretation of galaxy spectral energy distributi ons (for sim- 
ilar conclusions, from a modelling perspective, see iFerrara et alj 
1 19991 ; iGranato et ai]|2000l; iPanuzzo et alJl2007l;lRocha et alj | 2008l; 
from an observa tional perspective, see iBoquien et alj|2009l 120121 ; 
iMao et al.ll2012h . 



4.5 Summary 

In conclusion, therefore, the approach presented in this paper 
should help us gain valuable insight into the content and spatial 
distribution of stars and dust from spectral analyses of structurally 
unresolved star-forming galaxies. This is enabled by the possibil- 
ity to associate stars in different age ranges with different geometric 
components of sophisticated models of radiative transfer. It is worth 
emphasising that the results obtained in this way should not depend 
critically on the details of the adopted radiative transfer model, as 
the predictions of the versatile T04 model adopted in our MCMC 
analysis are similar to those of other types of state-of-the-art ra- 
diative transfer models (Section 12.21 and |2"3V We have illustrated 
the usefulness of this approach to interpret the dependence of in- 
tegrated spectral properties on galaxy inclination in a sample of 
nearby star-forming galaxies. Our approach also opens new possi- 
bilities to straightforwardly implement the effects of different spa- 
tial components and galaxy inclination in the predictions of semi- 
analytic models of galaxy evolution. This will be particularly useful 
to gauge the importance of selection effects and systematic biases in 
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the derivation of physical parameters from large samples of galax- 
ies using such models. 

We can summarise our main results as follows: 

(i) Four different types of radiative transfer models all predict 
a quasi-universal relation between slope of the attenuation 
curve at any wavelength from the ultraviolet to the near in- 
frared, nJ^ M (f?), and V-band attenuation optical depth of the 
dust in the diffuse ISM, t l y M (S), at all galaxy inclinations 
B (Section|2~2l. 

(ii) We parametrize the quasi-universal relation between 
»], SM (0) and f'y M (#) by mean of the simple analytic ex- 
pression given in equation ([9j of Section |272l This expres- 
sion can be easily implemented in any spectral analysis of 
structurally unresolved disc galaxies to account in a sim- 
ple yet physically consistent way for the effect of galaxy 
orientation on dust attenuation. 

(iii) The relation between rff u {6) and t l ^ M (8) given in e qua- 
tion ([9]( predicts steeper curves than the Calzetti ( 200 1|) and 
Milky- Way curves at small f'^ M (0), and shallower curves 
at larger f l y M (9). This may be rele vant to the in t erpret ation 
of recent observational results by Wu vts et"aT1 d201ll) and 
Arnouts et. al (in preparation), who find a systematic flat- 
tening of the optical attenuation curve with increased star 
formation activity in distant galaxies. 

(iv) We propose an original, computationally inexpensive ap- 
proach to combine the output of sophisticated radiative 
transfer models with models of spectral evolution, by as- 
sociating stars in different age ranges in the ambient ISM 
with different geometric components of galaxies (Sec- 
tion [3X2l . 

(v) We find that the comb ination of the radiative transfer 
model o flTuffs et alj ( 120040 with the galaxy spectral evolu- 
tion code of Chariot & Bruzual (in preparation) accounts 
remarkably well for the dependence of the Ha/H/3 ratio 
and ugrizYJH attenuation curve on axis ratio (i.e. inclina- 
tion) in the Wll sample of nearby star-forming galaxies 
(SectionEU. 

(vi) The central face-on 6-band optical depth favoured by this 
analysis is t b± « 1.8 ± 0.2 (corresponding to an angle- 
average (f^ M ) e « 0.3), slig htly larger than the range 
0.5 < t b± < 1.1 derived by IXilouris etail d 19991) from 
the study of nearby spatially resolved galaxies , and smaller 
than th e constraint t b± =s 3.7 ± 0.7 inferred bv lDriver et al.l 
d2007Q from the analysis of a sample of spatially unre- 
solved nearby galaxies. 

(vii) The roughly equal central face-on 6-band optical depths 
of low- and high-/j, galaxies, which have similar size dis- 
tributions in the Wll sample, are consistent with the an- 
ticorrelation between dust-to-stellar mass ratio and stellar 
surface mass density inferred from the ultraviolet, optical 
and infra red spectral analysis of low-redshift star-forming 
galaxies JCortese et al]|2012l) . 

(viii) For the galaxies in the Wll sample, over 80 percent of 
the attenuation is characteristic of that affecting thin-disc 
stars in the T04 model, while the remaining 20 percent 
is more similar to that affecting thick-disc and bulge stars 
(Section[3~4t. 

(ix) The median unattenuated V-band luminosity ratio of thick- 
disc to thin-disc stars inferred from our analysis of struc- 
turally unresolved, randomly oriented galaxies in the Wll 
sample is in remarkable agreement with that derived from 



the analysis of a sample of spatially resolved , nearly edge- 
on galaxies bv lYoachim & Dalcan ton (2006). 
(x) The influence of galaxy orientation and of the spatial dis- 
tribution of dust relative to stars on the shape of the attenu- 
ation curve can have a large impact on broadband colours 
of galaxies, reaching up to 0.3-0.4 mag at optical wave- 
lengths and 0. 1 mag at near-infrared ones. 
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APPENDIX A: GAUSSIAN RANDOM PROCESS 
EXTRAPOLATION 

Al Method 

We describe here the algorithm adopted to infer the attenuation 
curves of the T04 thick disc and bulge models in the range 0. 1 < 
A < 0.45 pin from the region 0.45 < A < 2.1 [im, employing 
the T04 thin disc model, which is computed in the whole range 
0. 1 < A < 2.1 pm. To make such an educated guess, we appeal 
to a Gaussian Random Proces (e.g. iRasmussen & WilliamsluOOq ; 
Wandelt 20120 This is a generalisation of the Gaussian probabil- 
ity distribution, which determines the relative likelihood of random 



variables, to functions. A Gaussian Random Process hence mea- 
sures the relative likelihood of different functions, that is of differ- 
ent relationships among random variables. To ease the understand- 
ing of a Gaussian Random Process, we recall here the definition of 
function: a function is a relation between a set of inputs and a set of 
admissible outputs with the property that each input is related to ex- 
actly a single output. Analytic functions (e.g. polynomials, power 
laws, splines) are a particular class of functions commonly used 
to solve interpolation and extrapolation problems, while a Gaus- 
sian Random Process is a powerful non-parametric method able to 
learn the (non-linear) relationships among a finite set of inputs and 
outputs. In particular, by appealing to a Gaussian Random Process 
we can solve any problem of interpolation and extrapolation (i.e. of 
inference) at arbitrary points from a finite set of known points (e.g. 
observations, simulations), within a fully Bayesian framework. A 
Gaussian Random Process is fully determined by a prior (over func- 
tions), and by a covariance matrix, which encodes the relationships 
among the known points. The result of such an approach is a pos- 
terior distribution over functions, which naturally accounts for the 
uncertainties due to a finite sample and to different possible rela- 
tionships among the known variables. 

To apply the algorithm, we compute a library of n = 2500 
attenuation curves from the T04 thin disc model (i.e. the 'train- 
ing library'), shown in Fig. lAlfa . We consider 50 values of t b± 
and 6 uniformly distributed in the ranges 0.0 ^ t b± ^ 8.0 and 
^ cos 8 ^ 1, respectively, and 100 wavelength bins uniformly 
distributed in log A in the range 0.1 < A < 2.1 pm. We refer to 
any thin disc attenuation curve in the training library in the range 
0.1 < A < 0.45 and 0.45 < A < 2.1 ^m as X, and Yj respec- 
tively. Any Xj curve contains n x = 49 wavelength bins, and any Y, 
curve n Y = 5 1 bins. We then calculate the mean attenuation curves 
fix = and \i Y = Yi- Next, we compute the n Y x n Y 

sample covariance matrix CvJ'l 

c ? y = jE(y/-/4)(^ -/4)> (ad 

where the indices k, j = l,ny span the wavelength bins in the range 
0.45 < A < 2. 1 pm. Then, we calculate the n x x n Y sample covari- 
ance matrix Cxr 

4V = a4) Of (A2) 

where the indices I = i,n x and j = \,n Y span the wavelength bins 
in the range 0.1 < A < 0.45 and 0.45 < A < 2.1 pm respectively. 
Finally, to infer the attenuation curve for any T04 thick disc and 
bulge model in the range 0.1 < A < 0.45 pm, we calculate the 
attenuation curve Y* in the range 0.45 < A < 2. 1 pm, and obtain 
the posterior mean over functions Hx*\ Y , which we call X*, in the 
range 0.1 < A < 0.45 pm 

X*= t i xv = l x x + C xr C^(r-Hr). (A3) 

The curve X* is the most probable estimate of any thin disc and 
bulge curve in the range 0.1 < A < 0.45 pm, given the curve in the 
range 0.45 < A < 2. 1 pm and the T04 thin disc model. The poste- 
rior covariance, which is a detailed quantification of the uncertainty 
in this prediction method, can be computed as 

Cx*|F = Cxx ~ Cx Y C yy Cyx ■ (A4) 



10 The Rasmussen & Williams 12006) book is freely available online at 11 We adopt here the maximum likelihood estimate of the sample covari- 

http: //www.Gaussianprocess .org/gpml/chapters/RW.pdf. ance matrix, so the denominator is n and not n - 1. 
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In the next two sections we test the accuracy of this approach 
by computing a set of attenuation curves in the range 0.1 < A < 2. 1 
\im, then comparing the curves predicted in the range 0.1 < A < 
0.45 |J.m with the 'true' ones. 

A2 Tests of the accuracy of the method 

A2.1 Test on the T04 thin disc model 

For the first test, we compute 20 attenuation curves from the T04 
thin disc model in the range 0.1 < A < 2. 1 \im by randomly draw- 
ing 20 pairs of 6 and t b± in the ranges ^ cos 6 ^ 1 and 0.0 ^ 
t~b±^; 8, respectively. We then use the training library described in 
Section [ATI to extrapolate the curves in the range 0.1 ^ A ^ 0.45 
\im from the region 0.45 < A ^ 2.1 \im by means of equa- 
tion dA3t . Fig. lAlb shows a remarkable agreement between the 
inferred and 'true' curves, though the prediction becomes less ac- 
curate for curves very close to zero, where there is little information 
to constrain the extrapolation. 

A2.2 Test on the \Gordon et aZI \200j) model 

This first test restricts the accuracy of the method in that we use 
the same model to generate the test sample and to build the train- 
ing set. Hence, to better determine the accuracy of the method, 
we use a more flexible model, which allows us to generate test 
curves for model parameters exploring a larger range than the train- 
ing set. We build a library of extinction curves in th e whole range 
0.1 < A < 2.1 |im following the parametrization of iGordon et al.l 
(2003), and, as in the previous test, try to infer the curves in the 
range 0. 1 < A < 0.45 p.m from the region 0.45 < A < 2. 1 ixm, 
eventually comparin g the predicted curves with the 'true' ones. 
IGordon et alj d2003h consider a comprehensive set of measure- 
ments of ultraviolet-to-near infrared extinction curves along dif- 
ferent line of sights of the Large and Small Magellanic Clouds, 
adopting the following 6 parameters function to fit the observed 
curves 

1.086"?,, = (AT, + K 2 x + K 3 D(x, y, x a ) + K 4 F(x))E(B - V) 

+ R V E(B-V), (A5) 
where x = A' 1 , and 



F(x) = 0.5392(x - 5.9) 2 + 0.05644(x - 5.9) 3 . 

We use equation dA5l > to compute a library of 1000 extinction 
curves (i.e. the 'training library'), shown in Fig. IA2t . sampled uni- 
formly in log A in the range 0.1 < A < 2.1 p.m. We consider four 
sets of parameters corresponding to the average 'SMC Bar Sam- 
ple', 'SMC Wing Sample', 'LMC/LMC2 S upershell Sample' an d 
'LMC Average Sample' given in table 3 of IGordon et al.l J2003h . 
To accomplish this, we randomly select one of the four sets of pa- 
rameters and randomly draw each parameter from a Gaussian dis- 
tribution centred on th e parameter with v ariance equal to the error 
quoted in table 3 of Gord on et al.1 d2003l) . Next, as for the library 
of T04 thin disc curves, we compute the mean curves /u x and n Y 
in the ranges 0.1 < A < 0.45 and 0.45 < A < 2.1 respectively, 
and the sample covariance matrices C YY and C X y as defined in 
equations dAlt and dA2t . respectively. After, we compute 20 test 
curves by randomly drawing the parameters K\, K2, K3, K4, xo,y 



from Gaussian distributi ons with variance equ al three times the er- 
ror quoted in table 3 of IGordon et al. I d2003l) . Finally, we calcu- 
late the curves in the region 0.1 < A < 0.45 urn from the range 
0.45 < A < 2. 1 \im through equation dA3t . We highlight that in 
this case the test curves span a wider range in the 6-dimensional 
parameter space than the library used to infer the curves, since the 
parameters A", , Ki, K 3 , A4, xo, y can reach more extreme values be- 
ing drawn from Gaussian distributions with a larger variance than 
the one adopted to compute the training library. Fig. lA2b shows the 
results of the Gaussian Random Process inference: it is noteworthy 
how well the curves in the range 0.1 < A < 0.45 |J.m are predicted 
from the region 0.45 < A < 2.1 urn, though the accuracy decreases 
at A < 0.25 u.m. This has no impact on our analysis, since the bluer 
filter used in this work is the SDSS u, whose transmission drops to 
zero at 0.3 |J.m. Note also that predictions obtained with a single 
Gaussian Random Process are accurate for both curves with and 
without the ultraviolet bump at 2175 A. 

A3 Application of the method to the T04 thick disc and bulge 
models 

The previous tests demonstrate that the extrapolation method de- 
scribed above, based on a Gaussian Random Process, is a reli- 
able way to infer attenuation and extinction curves in the range 
0.1 < A < 0.45 |j.m from the region 0.45 < A < 2.1 \im, exploiting a 
training sample and the correlation existing among different wave- 
length bins. We therefore apply the algorithm to the set of thick disc 
and bulge curves in the ranges ^ cos 6 ^ 1 and 0.0 ^ r B± ^ 8 
computed by T04. We show in Figs. |A3fa ,b the curves predicted in 
the range 0.1 < A < 0.45 \im from the region 0.45 < A < 2.1 ixm 
by means of equation dA3t . 

APPENDIX B: AN IMPORTANCE SAMPLING-BASED 
METHOD TO TREAT SMOOTH SYSTEMATICS 

Bl Methods: Monte Carlo and Importance Sampling 
integration 

Any method to measure dust attenuation based on the ratio of the 
spectral energy distributions of two objects (e.g. stars, galaxies) as- 
sumes that the intrinsic spectral properties (i.e. un-attenuated) are 
equal. In this work, we measure the effect of dust from the ratio 
of the (normalised) mean spectral energy distribution of galaxies 
in different axis ratio bins to the face-on one. This traces dust at- 
tenuation if the mean stellar populations properties of the galaxies 
are constant across the axis ratio bins, since this guarantees that the 
un-attenuated spectral energy distributions are also constant. We 
check this hypothesis by comparing the bin-averaged and sample- 
averaged values of 12 + log(0/H), ^s, Z and Rgo, finding that they 
systematically vary as a function of the axis ratio, hence biasing 
the measurement of dust attenuation. To correct the biases, we de- 
velop an original method which uses the concepts and formalism 
of importance sampling. While in usual computations of average 
quantities each object has the same weight, this method allows us 
to assign a different weight to each object in each axis ratio bin 
accordingly to how likely it is to be drawn from a distribution com- 
mon to all the bins. In this way we are able to substantially decrease 
the impact of biases on the attenuation measurements. We recall 
here the basic concepts of Monte Carlo and importance sampling 
(a Monte Carlo technique) in the general case of ^-dimensional dis- 
tributions. 
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Figure Al. Test of the accuracy of the Gaussian Random Process method (section [aJ by means of the T04 thin disc model, (a) 'Training' library of 2,500 
T04 thin disc attenuation curve. Solid points indicate the mean curve, along with the 100 log-spaced wavelength bins, (b) Test curves from the T04 thin disc 
model. Solid lines indicate the 'true' curves in the whole range 0. 1 < A < 2.1 u.m, dashed red lines the curves inferred in the range 0.1 < A < 0.45 |j.m from 
the range 0.45 < A < 2.1 p through equation iA3) . In both panels, the vertical dashed line separates the ranges 0.1 < A < 0.45 and 0.45 < A < 2A [im. 
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Figure A2. Test of the accuracy of the Gaussian Random Process meth od (section [A) by means of the lGordon eTaT] J2003l) model, (a) 'Training' library of 
1,000 extinction curves from the parametrization of lGordon et al.lf2003l) . The solid points indicate the mean curve, along with the 100 log-spaced wavelength 
bins, (b) Test curves. Solid lines indicate the 'true' curves in the whole range 0.1 < A < 2. 1 |xm, dashed red lines the curves inferred in the range 0.1 < A < 0.45 
u.m from the range 0.45 < A < 2.1 u.m through equation jA3) . In both panels, the vertical dashed line separates the ranges 0.1 < A < 0.45 and 0.45 < A < 2.1 
u.m. 



We start by recall ing basic Monte Carlo integration, following 
the nomenclature of iRobert & Casellal d2004l) . Suppose we want to 
evaluate the integral 

E[h(X)] = \ dxh(x)f(x), (Bl) 

over the domain X, where h(x) is a function, f(x) a probability dis- 
tribution, E indicates the expected (i.e. mean) value, and<Y 6 R". In 
practical problems, it is often impossible to solve the above integral 
analytically, and numerical techniques which are effective in low 
dimensions (i.e. quadrature rules) become inefficient in highly di- 
mensional problems (the so-called 'curse of dimensionality'). The 
Monte Carlo method allows the estimation of the integral l IB 1 1 
through the simple expression 

1 N 

h N = -Yjh( Xl ), (B2) 



where h N is the Monte Carlo estimate of E[/i(X)] and N the number 
of random variates x -, drawn from the distribution f(x). The power 
of this method is that it transforms the operation of integration into 
summation. To apply equation (B2i , we only need to draw random 
numbers from the multi-variate distribution fix), then, by means of 
the Strong Law of Large Numbers, the above expression converges 
to the 'true' value E[/»(X)], with variance 

1 N i 
var(h N ) = J^J] [h(xd - h N f ■ (B3) 
(=i 

Fig. IB 1 1 shows an example of (uni-variate) Monte Carlo integra- 
tion for h(x) = x with Gaussian probability distribution f(x) = 
N(fj,cr 2 ), lei, and fi = 5.743, cr 2 = 0.248. In this case one can 
easily solve the integral ( IB 1 1 analytically, obtaining the solution 
E[X] = 5.743. To compute the Monte Carlo estimate of the integral 
E(jc), we draw N random numbers from the probability distribu- 
tion f(x) = AC(5.743, 0.248), then compute ~x N with equation JB2I ). 
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Figure A3. Results of the extr apolation o f the T04 thick disc and bulge models by means of a Gaussian Random Process, («) Solid lines indicate the T04 
thick disc model computed by Tuffs et al. 12004) in the range 0.45 < A < 2 um, dashed lines the curves inferred in the range 0.1 < A < 0.45 |^m through 
equation iA3\ . (b) Same as (a), but for the T04 bulge model. 



. E|xf = 5.743 ' 

■ X = 5.819 (N= 100) 

■ X = 5.709 (N = 1000) 

■ X = 5.749 (JV = 10000) 
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Figure Bl. Example of Monte Carlo integration of the univariate Gaussian 
function f(x) = Af(5.743, 0.248) (black line). Red, blue and dark green lines 
are histograms of 100, 1000 and 10000 samples drawn from the Gaussian 
distribution Af(5.743, 0.248), used to compute the integral of f(x) by means 
of equation <B2t . 



We show in Fig. IB II the random samples for N = 100, 1000 and 
10000, along with the corresponding x N . Note that the error on x N 
decreases approximately as 1 / VJV, as expected from equation < IB3b . 

In the previous example, we have solved the integral l lB 1 1 
in the simple case of x following a Gaussian probability 
distribution/fx) = N{p.,a 2 ), exploiting the fact that we can effi- 
ciently draw samples from a Gaus sian distribution, us ing, for in- 
stance, the Box-Muller transform dBox & Mullenll958l) . In many 
situations it is not feasible, or computationally efficient, to draw 
samples from f(x), so we can exploit the fact that the integral ABU 
is equivalent to 



E[h(X)] 



L 



dxh(x) 



fix) 



(x), 



(B4) 



where we have introduced the 'importance function' g(x), to 



rewrite the Monte Carlo estimate (IB2) as 



h 



Zm /(X.)/g(Xi) 



(B5) 



Equation dB5b expresses the importance sampling estimate of the 
integral JBU . and it differs from the Monte Carlo estimate l lB2t 
since it allows the evaluation of properties of a distribution (e.g. the 
expected value) by drawing samples from a different distribution. 
Also, importance sampling is a variance reduction technique, since, 
while equation JB2b implies a uniform sampling from the probabil- 
ity distribution f(x), equation dB5t allows one to concentrate the 
sampling on the most 'important' regions of f(x) (i.e. those that 
contribute the most to the integral), by drawing samples from a dif- 
ferent distribution g(x). Equation l |B5l > is similar to dB 21 >. except 
for the weight factor w(Xj) = f(xi)/g(Xi), which accounts for the 
fact that we draw samples from the 'wrong' distribution. Here we 
normalise the weights to 1, so the functions f(x) need not to be nor- 
malised. In this case the quantity h 1 ^ is biased with bias jGewekd 
1989) 



775-, 2m *(*/)/(*!)/*(*() 



and variance 



var(h<S) : 



Sh f(x,)/s(x,) 



N Zm(^(x/)-M 2 
JV- 1 



(B6) 



(B7) 



To be a valid importance function, g(x) has to fulfil the following 
condition 



g(x) = 



/(x) = 0, 



(B8) 



which guarantees the equality of the integrals < IB H and l !B4t . This 
means that the function g(x) has to be non-null for every x in the 
support of f(x), or, in other words, that all the values of x which 
have non-null probability in f(x) must have non-null probability 
also in g(x). Though this is a sufficient condition for an importance 
function, in practical applications we need a further constraint, that 
is an estimator with finite variance. This translates in a further con- 
dition 



for M > and Vx, < M 
8(x) 



var(h N ) < oo , 



(B9) 



© 0000 RAS, MNRAS 000, 000-000 



Content and spatial distribution of dust in galaxies 29 



— E[X]'= 5.743 ' 

— X = 5.803 (N = 100) 

— X = 5.736 (N = 1000) 

— X = 5.743 (JV = 10000) 




Figure B2. Example of Importance Sampling integration of the univariate 
function f(x) = N (5.743, 0.498)0.2 sin x 2 + JV(5.743, 0.498) (black line). 
Dot-dashed line is the importance function g(x) = Af(5.743, 0.77), red, blue 
and dark green lines are histograms of 100, 1000 and 10000 samples drawn 
from the importance function g(x) used to compute the integral of f{x) by 
means of equation IB41 . 



which means that, for every x in the domain X, the tails of g(x) 
must be at least as heavy as the ones of f(x). Fig. |B2] exempli- 
fies the importance sampling integration of the uni-variate function 
f{x) = W(5.743, 0.498) 0.2 sin x 2 + W(5.743, 0.498). Since we can 
not easily draw samples from f(x), we adopt a Gaussian distribu- 
tion g(x) = N(p. x , oy) as importance function and draw samples 
from it. We then compute the importance sampling estimate of f(x) 
with equation JB5t . obtaining the estimates reported in Fig. lB2l 



B2 Application of the method to the biases of the Wll 
sample 

We now describe how we make use of the importance sampling 
formalism to correct for the biases in the dust attenuation measure- 
ments. We can arrange the gas-phase metallicity 12 + log(0/H), 
specific star formation rate t/^, redshift z, and radius which con- 
tains 90 per cent of the r-band Petrosian flux i?9o, in a vector 
x = [12 + log(0/H), tf/ s ,z,Rgo]- Each galaxy is then associated to 
a vector Xj, and we can write the systematic variation of the bin- 
averaged values of 12 + log(0/H), tf/s, z and Rw with the axis ratio 
b/a as 

E[p(x\b/a)] = E[p(x)] + A(x\b/a), 

where A(x\b/a) indicates the axis ratio dependent bias. We can 
re-write the above expression by using the definition of expected 
value jB I | i in the case of a discrete distribution 

Xip[xi\(b/a)i] = ^ Xip(xi) + A[Xi\(b/a)i] , 

i i 

where we sum over al possible values of x, and p{xt) = p(xj) if 
we assume that the observations are independent and equiproba- 
ble when accounting separately for the effect of the axis ratio de- 
pendent bias. We can now introduce a new probability distribution 
p'(x), which we assume independent of the axis ratio b/a, so as 
p'(x\b/a) = p'(x). We can interpret this as the probability distribu- 
tion f(x) of equation JB5I >. and p(x \b/a) as the importance function 



g(x) of the same equation dB5t . If we apply equation JB4b . replac- 
ing f(x) with p'(x\b/a) and g(x) with p(x\b/a), we obtain 

-_ 'LiXip'{x i )lp[x i \(bla) i \ 
X ~ ZiP'OrO/ptel^/a);] ' 
which is by construction constant across the axis ratio bins. 

To apply this method we must be able to compute the proba- 
bility distribution p(x\(b/a)) in each axis ratio bin. To accomplish 
this, we can not rely on multi-dimensional histograms, since the 
probability distributions so obtained are non-smooth and depen- 
dent on the choice of the bin widths, hence we appeal to a kernel 
density estimation. This uses a kernel to assign to each point in the 
space defined by the parameters x a density (i.e. a probability). A 
commonly adopted kernel is the Gaussian kernel, which weights 
each data point with a multi-variate Gaussian function centred on 
the point itself, with diagonal covariance. This allows a smooth es- 
timate of the joint probability distribution, though the choice of the 
diagonal covariance is critical. Fortunately, it is possible to find an 
optimal covariance, also called filter bandwidth, through the mean 
square error minimisation, by adapting the Gaussian kernel covari- 
ance to the density of the points at each location ('a daptive Gaus- 
sian kernel dens ity estimato r', see ISilverman|[l986T) . We use the 
LibAGF library (Mills 2011) to compute for each axis ratio bin k 
a smooth probability distribution p^(x\(b/a)k), which can be evalu- 
ated at arbitrary points in the parameter space defined by the vector 
x. Next, we need to find a suitable probability distribution p'(x), 
so as criteria ( |B8) dB9t are valid. We therefore adopt the function 

p'(x) = Y]pk(x\(b/a) k ) 0A , 

which satisfies condition dB8t . since a null value in any of the func- 
tions p k (x\(b/a) k ) make null the corresponding value of p'(x). The 
choice of the exponent of the product of the densities p k (x\(b/a) k ) 
determines, for each axis ratio bin, the fraction of galaxies with 
significant weights; an exponent equal to 1 produces a very peaked 
probability distribution p'(x), so as less than 1 % of the objects 
within each bin have significant weights. This produces an almost 
perfect erase of the biases, assuming that the physical quantities 
x = [12 + log(0/H), i//s,Z,Rgo] are perfectly measured. Otherwise, 
this result would be dominated by biases and errors on the measure- 
ments of the physical parameters x. After some experimentation, 
we find that a value of the exponent of 0. 1 produces a decrease of 
the biases at the level of < 10 percent, while keeping the fraction 
of objects with significant weights above 10 percent, and hence de- 
creasing the risk of being dominated by measurement biases and 
errors. To check the validity of condition ( TB9) , we examine the be- 
haviour of the tails of the distribution p'(x ) through the histograms 
of the weights w k (Xj) = p'(x)/p[Xi\(b/a) k ]. In the end, we obtain 
for each galaxy in the low- and high-yu, samples a weight w k (Xi), 
which we use to compute the expectations (i.e. the mean values) of 
the Ha/H/? ratio and ugrizYJH attenuation in each axis ratio bin. 



APPENDIX C: MARKOV CHAIN MONTE CARLO 

In a Bayesian framework, the degree of belief in a set of parameters, 
which define a model, given a set of data, is quantified by the poste- 
rior probability distribution defined by Bayes' Theorem. Evaluating 
the posterior probability distribution from a multi-dimensional grid 
of parameters can be computationally demanding because, at fixed 
accuracy, computational time scales exponentially with the number 
of parameters. To overcome this limitation, the posterior probabil- 
ity distribution can be represented through a set of samples drawn 
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from it. In practice, a sampled representation allows an easy com- 
putation of one-dimensional (or two-dimensional) marginal pos- 
terior probabilities, which provide a quantitative summary of the 
knowledge about individual parameters (or different pairs of pa- 
rameters) given a set of data. 

Markov Chain Monte Carlo (MCMC) is a general method to 
generate random samples from a large class of probability distri- 
butions for which direct sampling alg orithms are not k nown. Sev- 
eral MCMC algorithms exist (see iGelman et alj|2010l . and refer- 
ences therein), some based on a random walk through the param- 
eter space. The power of random walk-based MCMC algorithms 
is that the set of all points visited by the random walk converges 
towards the desired target distribution. In this work, we consider as 
the target distribution the posterior probability distributions of the 
five adjustable parameter of the dust model t^ c ,T BJ _ ,t t hi n , t t hkk and 
tbuige (see Section [3.2.2^ In the following section we will give de- 
tails of the specific random walk method adopted, and discuss the 
diagnostics that we use to asses the convergence of the MCMC to 
the target distribution. 



CI The Metropolis-Hastings algorithm 

The M etropolis-Hastings sampler dMetropolis et al .ll953l : lHastingj 
Il97fj|) is one of the most widely used MCMC algorithms. The algo- 
rithm iteratively explores the parameter space by means of a ran- 
dom walk, eventually converging towards the target distribution. 
It needs a set of starting parameters Xj =0 , which are usually ran- 
domly drawn within the priors range. Then, at each iteration j, a 
new set of parameters X', is drawn from a proposal distribution 
q(X'j\X and the relative probabilities of the parameters X'j and 
Xj-i are compared by means of the Metropolis ratio, defined as 

PiX'p,!) g(X hl \X } ) 
p(Xj-,\D,I) qiX'jlXj-i) 

Next, the Metropolis ratio r is compared with a random number n 
uniformly drawn between and 1, and the new set of parameters 
is accepted if n ^ r (setting Xj = X'j), otherwise is rejected (set- 
ting Xj = Xj-i). In the case of a symmetric proposal distribution, 
the second factor of equation ( IClt is 1, therefore the Metropolis 
ratio reduces to the likelihood ratio for the two set of parameters 
X'j and Xj-\. To obtain a stationary target distribution, these steps 
must be repeated for a long enough time, which depends on the 
way the parameter space is explored. This exploration is set by the 
adopted proposal distribution, which determines the moves in the 
parameter space. The proposal distribution has to satisfy the con- 
dition that any point in the support of the target distribution can be 
reached from any other point through a set of proposal moves. The 
commonly adopted Gaussian proposal distribution fulfils this con- 
dition, since it has infinite support, though in the case of pathologi- 
cal posteriors it can fail in the exploration of some regions of the pa- 
rameter space. The optimal proposal distribution, which minimizes 
the number of moves and the computational time, has the same 
shape of the posterior distribution. Though, in many situations this 
information has little relevance, since the posterior distribution is 
not known a priori, being the goal of the algorithm to compute 
it. A solution is to use an adaptive proposal distribution, which 
has the property of learning the shape of the posterior distribution 
as the chain moves through the parameter space (see Haa rio et alj 



l200lHBrooksetalJl201lL and references therein)^ The package 
CosmoMC implements an adaptive Gaussian proposal distribution, 
whose covariance is computed from the posterior samples obtained 
so far. In our work we have no previous knowledge of the posterior 
probability distribution, hence such an adaptive algorithm needs to 
(inefficiently) obtain a lot of samples before obtaining a good es- 
timate of the posterior covariance, which is in turn used to set the 
covariance of the Gaussian prop osal distribution. H ence, we use 
the nested sampling algorithm (Skirling 2004, 2006) to rapidly ex- 
plore the parameter space and sample the unknown posterior dis- 
tribution, then implementing in CosmoMC the posterior covariance 
obtained with such an algorithm. Also, this allows us to check the 
consistency of the posterior distribution obtained with two different 
algorithms. 

C2 Convergence diagnostics of the MCMC run 

The main difficulty of the Metropolis-Hastings algorithm, and in 
general of random-walk based MCMC, is the determin ation of the 
convergence of the chain to the target distribution (see iKass et al] 
ll998l ; ICowles & Carlinl 19961 and references therein). A way to ac- 
complish this is to run multiple parallel chains, comparing on-the- 
fly the samples obtained from the different chains. CosmoMC im- 
plements such a convergence algorithm, comparing the between- 
and within-chain variances from the last half of the sam ples of 
each parameter, by means of the IGelman & Rubin! |l992) R pa- 
rameter. In this work, we run in = 6 chains, hence after n steps 
we have, for each parameter, n x in samples, that we indicate as 
Xj h j = 1, ...,m; i = 1, To compute R, we consider, for each 
parameter X separately, the between-chain variance B defined as 

where 

1=1 .;=! 

are the mean of the parameter X within any chain j and across all 
the chains, respectively. We then define the within-chains variance 
Was 

1 "' 

w=-y S] 2 . 

j=i 

where 

is the variance of the parameter X in the chain / Finally, we can 
compute the convergence parameter R as 

In- 1 TIT 

R= J + -— . 

V n nW 

A value of R < 1 . 1 for each parameter of interest usu ally indicates 
good convergence (see Section 1 1.6 of lGelmanll2004h . In our case, 
we set in CosmoMC the converge criterion R ^ 1.01. 

To further check the exploration of the parameter space, we 

12 Chapter 4 of [Brooks et al] | |201 lb is also freely available online 
at http : // citeseerx . ist . psu . edu/viewdoc/download?doi=18 . 1 . 
1 . 161 . 2424&rep=repl&type=pdf 
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plot in Figs. IClb4C2b the value of the draw of each parameter at 
each iteration against the iteration number (i.e. the so called 'trace- 
plot') for three randomly selected chains among the six chains run. 
Figs. lClfa4C2fa show that both for the low- and high-/j„ data sets the 
chains move from the starting points to the region of the parameter 
space of highest probability, without showing any large-scale trend 
which could suggest a bad exploration of the parameter space. This 
also indicates that discarding the first half of any chain, as sug- 
gested in CosmoMC, is a conservative choice which allows us to 
keep only draws that are likely sampling the target distribution. 

Finally, we plot in Figs. IClb4C2b the value of the auto- 
correlation function between any pair of draws as a function of 
the separation between the draws (i.e. the 'lag'), for each parame- 
ter and for the same three randomly selected chains shown in Figs. 
|Clb - IC2b . The auto-correlation quantifies the amount of correlation 
between subsequent draws, and starts from 1 at lag 1 = 1, since each 
move depends on the previous one in a deterministic way through 
the proposal distribution. In the case of a good exploration of the 
parameter space, the auto-correlation function drops to ss after a 
lag /', meaning that samples drawn at iterations ( and i + I' are in- 
dependent. Figs. |Clb4C2b show that the auto-correlation function 
drops and starts oscillating around zero at /' ss 70 for both the low- 
and high-//, data sets, as expected in the case of good exploration 
of the parameter space. 

These tests indicate that the chains are well exploring the pa- 
rameter space, and that the drawn samples have likely converged to 
the stationary target distribution, which in this case is the posterior 
probability distribution of the five parameters t^ c ,T B± ,t thin , t thick and 

tbulge • 



ity distribution. In practice, this is accomplished by computing, at 
each iteration, the likelihood on N points uniformly drawn from 
the prior distribution tt(0), then selecting the point with the low- 
est likelihood £. min and, at the next iteration, drawing points from 
7r(0) conditioned to have _£(0) > -£,„;,,. The technical difficulty 
is to efficiently draw points from the conditional prior distribu- 
tion tt(0|.£(0) > -Cm,,), since the volume occupied by parame- 
ters for which £(0) > -£„„■„ shrinks also exponentially with the 
iterations. The MultiNest algorithm faces this problem by appeal - 
ing to the 'simu lt aneo us ellipsoidal nested sampling method' of 
iFeroz & Hobsonl ( 12008b . This allows the approximation of con- 
tours with constant likelihood with ellipsoids, whose optimal num- 
ber and bounds are determined by an expectation-maximisation 
algorithm. This allows an efficient sampling from the conditional 
density tt(0|X(0) > -£,„;„), hence solving the major technical ob- 
stacle of the nested sampling algorithm. Also, MultiNest exploits 
the possibility of parallelisation of the nested sampling algorithm, 
strongly reducing the wall-clock time required for inference on a 
mulit-core computer. 

In this work, following authors recommendations given in the 
manual of the package, we run MultiNest with N = 500 active 
points, setting the tolerance of the evidence, which is the stopping 
criterion, to 0.5, and the sampling efficiency, which determines the 
balance among the accuracy of the posterior and the length of the 
run, to 0.7. The algorithm converges after about 15,000 likelihood 
evaluations. We then compute the covariance among of all pairs of 
parameters, and include it in the CosmoMC run. As a further consis- 
tency check, we compare the marginal and joint posterior probabil- 
ity distributions obtained with MultiNest and CosmoMC finding 
an excellent agreement among the two. 



APPENDIX D: NESTED SAMPLING 

As discussed in Appendix [C] the Metropolis-Hastings algorithm 
with Gaussian proposal distribution has two major limitations: the 
tuning of the proposal distribution, which is crucial to guarantee 
a good exploration of the parameter space, and the assessment of 
the convergence of the chain, which is non-trivial in the general 
case of multi-modal distributions. These two caveats become even 
more serious during the very first exploration of a model against a 
set of data, since the posterior distribution is unknown and the pro- 
posal distribution can not be tuned in an optimum way. This usually 
makes the first runs of a MCMC computationally demanding, since 
many likelihood evaluations are needed before one can start tuning 
the proposal distributions using the samples drawn from the poste- 
rior. Moreover, a single MCMC can fail in the exploration of multi- 
modal posterior distributions with well separated peaks, so mul- 
tiple chains and/or annealing (or tempering) schemes are needed 
to ensure the exploration of all modes. Another limitation arises 
from the nature of Markov Chains (i.e. the markovian property), 
which makes each step in the chain dependent on the previous one, 
making impossible to parallelise single MCMC. Motivated by the 
need of a fast and reliable a lgorithm to compute the Bayesian ev- 
idence, Skillin3 l200l, l2006) proposed a novel Monte Carlo ap- 
proach, which allows the computation of the multi-dimensional ev- 
idence integral in Bayes's theorem by transforming it into a one- 
dimensional integral. The algorithm, called nested sampling, ex- 
plores the parameter space in a stochastic way by means of set of 
'active points', shrinking exponentially at each iteration the prior 
volume from which the samples are drawn around the regions of 
highest probability. Besides the evidence, the drawn samples al- 
low, as a by-product, the computation of the posterior probabil- 
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Figure CI. Exploration of the five-dimensional parameter space defined by fy C , tb±, khin, khick and tbuige by the MCMC described in Section l33l (a) Iteration 
number against value of the draw of each parameter, at each iteration (i.e. 'trace-plot') for three randomly selected chains among the six chains run (indicated 
by the different colours), for the 1ow-/j, galaxies. Vertical dashed lines separate the samples discarded (first half of each chain) from those used to compute the 
posterior probability distribution (second half of each chain), (b) Auto-correlation function as a function of the lag (i.e. the distance between two iterations) 
for each parameter, for the last half of the same three randomly selected chains of (a) (indicated by the different colours). 
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Figure C2. Exploration of the five-dimensional parameter space defined by f^ c , tb±, t t hin, t t hick and tbuige by the MCMC described in Section l33l (a) Iteration 
number against value of the draw of each parameter, at each iteration (i.e. 'trace-plot') for three randomly selected chains among the six chains run (indicated 
by the different colours), for the high-/j, galaxies. Vertical dashed lines separate the samples discarded (first half of each chain) from those used to compute 
the posterior probability distribution (second half of each chain), (b) Auto-correlation function as a function of the lag (i.e. the distance between two iterations) 
for each parameter, for the last half of the same three randomly selected chains of (a) (indicated by the different colours). 
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